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PREFACE 


The  first  volume  of  this  monograph  included  the  formu¬ 
lation  and  self-similar  solutions  to  blast  waves.  The 
next  step  in  our  exposition  of  Blast  Wave  Theory  and  Compu¬ 
tations  then  expresses  the  need  for  presenting  analytical 
solutions  that,  in  effect,  extend  che  range  of  validity  of 
self-similar  ones.  The  authors  found  that  the  one  analytical 
method  that  is  most  amenable  to  solving  different  classes 
of  problems  is  the  Quasi-Self-Similar  Method  developed  by 
Oshima.  It  was  then  decided  to  devote  the  second  volume 
of  this  series  to  presenting  a  complete  exposition  to  this 
method  that  includes  solutions  to  several  types  of  problems, 
most  are  included  here  for  the  first  time. 

N,  ' 

-~‘>Tn  this  volume,  after  a  concise  introduction  in  Chapter 

I,  the  basic  blast  wave  equations  are  presented  in  Chapter 

II.  A  novel  treatment  of  the  quasi-similar  formulation, 
presenting  it  as  a  zero  order  solution  to  a  double  Taylor's 
series  expansion  is  presented  in  Chapter  III.  Applications 
are  then  presented  in  Chapter  IV  which  is  divided  into  six 
parts.  The  first  two  parts  include  solutions  to  the  adiabatic 
point  explosion  problem,  the  latter  being  an  analytical  closeu 
form  solution.  The  third  part  deals  with  blast  waves  m 

real  gases,  where  the  formulation  is  presented  here  in  more 
details  than  that  presented  in  Volume  I.  Blast  waves  i~  a 
detonating  medium  is  included  in  the  fourth  part,  with  energy 
added  in  the  front.  Blast  waves  in  reactive  media  where  the 
detonation  energy  is  released  in  a  spacially  varying  fashion, 
are  given  in  the  fifth  part.  Finally,  in  che  sixth  part, 
the  effects  of  viscosity,  heat  conduction  and  radiation  on 
detonation  are  presented.  .. 

Ths  authors  are  greatly  indebted  to  Mrs.  Ellen  Kamel, 

Mrs.  Janie  Abdel  Aziz  and  Mr.  Taher  Hour  for  their  ability 
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CHAPTER  I 


NON-SELF-SIMILAR  BLAST  WAVES 


1.1  INTRODUCTION 

Blast  waves  considered  here  are  geometrically  symmetrical,  non-steady  -flow  field 
of  a  compressible  medium  that  are  bounded  by  gasdynamic  discontinuities.  Gene¬ 
rally  they  are  formed  by  explosions.  The  process  is  governed  by  spatially  one¬ 
dimensional  time  dependent  equations  expressing  the  conservations  uf  mass,  momen 
turn  and  energy,  subject  to  the  appropriate  boundary  conditions  at  the  center  and 
at  the  front  for  the  particular  problem  under  consideration. 

1.2  SELF-SIMILAR  BLAST  WAVES 

In  Volume  I  of  this  report,  it  has  been  shown  that  under  certain  conditions  the 
blast  wave  equations  may  be  reduced  to  a  set  of  ordinary  differential  equations 
in  terms  of  a  single  similarity  variable,  thus  making  them  amenable  to  simple 
analysis.  This  may  occur,  in  general,  when  the  dependence  of  the  gasdynamic  pa¬ 
rameters  on  the  front  coordinate  is  anihilated.  Such  self -similar'  problems  are, 
as  a  rule,  characterized  either  by  a  constant  front  velocity,  or  by  a  negligible 
essentially  zero,  counter-pressure  that  causes  the  Mach  number  of  the  wave  to 
remain  infinite,  irrespective  of  its  actual  velocity.  In  addition,  the  specific 
heat  ratio  of  the  gas,  as  well  as  any  existing  source  terms  of  mass,  momentum  or 
energy,  may  only  be  admitted  js  functions  of  the  non-dimensional  similarity  va¬ 
riable.  If  the  wave  front  were  moving  into  a  variable  density  atmosphere,  then 
the  ambient  density  dependence  should  be  restricted  only  to  a  constant  power  of 
the  distance  from  the  center  of  symmetry. 

The  power  and  versatility  of  the  concept  of  self-similarity  have  been  demonstra¬ 
ted  in  Volume  I  by  a  wide  variety  of  possible  applications. 

In  many  physical  situations,  however,  especially  when  one  is  interested  in  the 
interpretation  of  experimental  records,  wave  front  Mach  numbers  are  finite  ant! 
the  waves  either  accelerate  or  decay.  Under  sucn  circumstances  similarity  con¬ 
ditions  are  inapplicable.  One  then  has  to  take  into  account  the  dependence  of 
the  gasdynamic  parameters  cf  the  problem  on  the  change  in  conditions  at  the 


front,  as  well  as  their  change  within  the  flow  field.  It  thus  becomes  necessary 
to  contend  with  a  system  of  non-iinetr,  coupled,  non-homogeneous  set  of  partial 
differential  ec nations. 

1.3  THE  OPAJI-SIMILAR  ’4FTEOD 

An  analytical  roluc.nn  to  this  set  of  equations  utilizing  the  so-called  "quasi- 
similar"  method  developed  by  Oshima  fl96Q,  1?62,  1964).  By  postulating  a  'sepa¬ 
ration  of  variables"  relation  between  the  dependent  and  the  independent  varia¬ 
bles,  then,  due  to  the  logarithmic  nature  of  the  governing  equations,  all  the 
terms  containing  the  front  cooroinate  may  take  their  values  at  the  front.  The 
equations,  thus,  are  reduced  to  a  system  of  ordinary  differential  equations, 
identical  to  the  self-similar  ones,  with  the  dependence  on  the  front  coordinate 
included  only  as  an  additive  algebraic  term  that  has  the  front  coordinate  as  a 
parameter.  An  explanation  of  the  theoretical  basis  of  the  quasi-similarity  me¬ 
thod  was  given  later  by  Kamel  (1971).  Ha< proved  that  the  quasi-similarity  ap¬ 
proximation  is  the  zerotn-ordar  solution  of  a  double  Taylor's  series  expansion 
of  the  dependent  variables. 

The  numerical  solution  to  the  quasi-similar  equations,  for  constant  energy  blast 
waves,  has  been  tabulated  by  Lewis  (1961)  for  the  plane,  cylindrical,  and  sphe¬ 
rical  syrrvnetries  with  a  variety  of  specific  heat  ratios,  while  Oshima  (I960, 
1964)  obtained  an  approximate  closed  form  solution.  This  method  was  then  used 
by  Lee  (1S67)  to  study  the  non-uniform  propagation  of  imploding  shocks  and  deto¬ 
nations,  while  Oshima  (1S67)  attempted  to  apply  it  to  the  decoupled  shock-defla¬ 
gration  system.  Rae  (1963,  1965,  1968),  on  the  other  hand,  used  it  in  conjunc¬ 
tion  with  his  meteoroid  impact  solutions.  Attia  (1974)  used  this  method  to  ob¬ 
tain  a  non-self-similar  solution  of  blast  waves  in  adiabatic  point  explosions 
based  on  the  criterion  of  zero  particle  velocity  at  the  cental  of  symmetry.  In 
Chapter  XI  of  Volume  I  of  this  report,  the  quasi-similar  method  was  used  to  ob¬ 
tain  a  solution  to  blast  waves  taking  into  account  real  gas  effects  (Ghuneim, 
1975 j  Kamel  et  al.,  1977b).  The  problem  of  non-self-similar  detonation  waves 
was  solved  by  Abdel-Raouf  (1982)  and  that  with  transport  chenomena  taken  into 
account  was  also  solved  by  the  same  author  (1982)  Qhyagi  et  al.  (1981)  took 
into  account  the  distribution  of  fuel  concentration  in  a  reactive  medium  by  an 
exponential  model  of  heut  release,  in  terms  of  the  front  coordinate,  and  ob- 
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tained  a  solution  for  the  spherical  wave  by  utilizing  the  quasi-similar  tech¬ 
nique  . 

In  this  volume,  the  quasi-similar  method  of  Oshima  Is  analyzed  and  discussed 
briefly,  since  it  is  a  simple  and  qualitatively  accepted  analytical  method  of 
solution  which  has  a  wide  range  of  applicability.  A  wide  range  of  possible 
applications  are  included.  These  applications  deal  with  blast  waves:  In  in 
adiabatic  medium)  with  conduction  and  radiation  effects  and  with  real  gas  ef¬ 
fects,  as  well  as  flame  and  detonation  induced  blast  waves. 


CHAPTER  II 


BASIC  EQUATIONS 


II. 1.  INTRODUCTION 

The  phenomenon  of  wave  propagation  in  gases,  constitutes  unsteady  motion  of 
the  flowing  medium,  which  is  governed  by  the  conservation  equations,  namely 
the  conservation  of  mass,  momentum,  and  energy.  When  the  condition  of  sym¬ 
metry  is  satisfied,  the  flow  field  behind  the  shocK  wave  is  one-dimensional, 
i.e.,  the  gasdynamic  parameters  will  then  depend  only  on  one  spatial  coor¬ 
dinate  \T  ,  measured  from  the  center  of  explosion,  and  the  time  t  ,  mea¬ 
sured  from  tf id  instant  of  explosion. 

In  this  chapter,  the  fundamental  conservation  equations  are  formulated  appli¬ 
cable  for  the  three  geometrical  symmetries;  planar,  cylindrical  and  spheri¬ 
cal,  without  restrictions  imposed  by  a  specific  form  of  an  equation  of  state 

* 

and  including  any  possible  source  terms  of  mass,  momentum  or  energy  that  may 
affect  the  flow  field.  The  fundamental  equations  are  thus,  accordingly, 
transformed,  with  the  use  of  appropriate  non-dimensional  variables  and  para¬ 
meters,  to  a  most  concise  non-dimensional  form,  of  Eulerian  space  profiles, 
for  an  arbitrary  equation  of  state  relating  the  internal  energy  with  pressure 
and  oensity.  The  formulation  of  the  problem  is  then  completed  by  the  speci¬ 
fication  of  the  boundary  conditions  imposed  by  a  gasdynamic  discontinuity  at 
the  wave  front. 

The  integral  relations,  which  are  the  global  conservation  equations  of  mass, 
momentum  and  energy,  are  also  formulated  in  their  general  form.  These  rela¬ 
tions  serve  a  dual  purposej  in  addition  to  providing  a  enecK  on  the  accuracy 
of  the  flow  field  obtained  by  solving  the  conservation  equations,  they  aleo 
provide  means  for  the  determination  of  the  front  trajectory. 

Finally,  the  reduction  of  the  conservation  equations,  boundary  conditions  and 
integral  relations  to 'their  simplest  forms  for  a  medium  which  behaves  as  a 
perfect  gas  with  constant  specific  heats  as  well  as  In  she  case  of  strong  ex¬ 
plosion,  or  self-similar  motion,  is  determined. 
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II. 2.  THE  CONSERVATION  EQUATIONS 

The  governing  equations,  when  taking  into  account  possible  sources  of  mass, 
momentum  or  energy  that  may  affect  the  flow  field,  are  expressed,  as  given 
by  Oppenheim  et  al.  (1971)  as  follows: 

» 

(II. 


P 

r 


vr 


SL 


Lr> 


(II. 2) 


[<*♦  Jj?  )  tr*  j  +  pf  +  £  +  P  )  tu  rJ]  *  /rt/2.  (n. 


3) 


where  U  is  the  particle  velocity,  f  is  tha  density,  P  is  the  pressure  and 
6  is  the  specific  internal  energy.  The  index  J  is  a  geometrical  factor  de¬ 
fined  in  terms  of  the  variation  of  the  flow  cross-sectional  area  A  with  di¬ 
stance  as: 


im  JJli.  =  O  J  i,  2 

j/dr 

for  plane,  cylindrical  and  spherical  symmetrical  flow,  respectively,  rhe  sym- 
dcIs  yi ^  ,  Sip  and  Sl^,  represent  the  mass,  momentum  and  energy  suurcj  terms, 
respectively,  per  unit  mass  of  the  flowing  medium. 


11.3.  THE  EQUATION  OF  STATE 
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Th8  governing  equations,  Eqs.  (II. 1)  -  (II. 3),  are  completed  by  the  equation 
of  state  that  expresses  the  internal  energy  in  terms  of  pressure  and  density. 

Thus,  since  €  »  C  (  P  ,  f  )  it  follows  that 

dl d  AiP  ^($5-5  )  (II. 4) 

from  the  above  equation,  one  has 


where  subscript  S  denotes  entropy. 

While  for  a  constant  Bntropy  process,  the  first  law  of  thermodynamics  requi 
res  that 


-  dP 

“  ~r 


where  L  is  the  specific  enthalpy. 
Thus 


or 


Je  -  PJf 
S’1 


from  the  above  equation,  one  has  in  general 
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(  ~og  \  _  _P  /  \ 

Hence  it  follows  that 


defining  now  a  non-dimensional  isentropic  compressibility,  or 
sound  modulus 


r 


\  _  9  a} 

C"^Tp  s  "  P 


and  a  non-dimensional  internal  energy  factor 


~og  \ 


Thus,  Eq.  (II. 5)  reduces  to 


<  P 

r  ' 


and  Eq,  (I I. 4)  becomes 


J«  «  JU-t  |  itnf 

r  t  f 


(II. 5) 


velocity  of 


(II. 6) 


:  1 1 . 7 ) 


(ii. a) 


til.  a 
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Equation  (II. 9)  represents  the  general  equation  of  state  in  a  differential 
form. 

Taking  partial  derivatives  of  the  above  with  respect  to  the  space  coordinate, 
IT  ,  and  the  time,  t  .  respectively,  one  obtains  finally: 


ne  —  X  +  i  P  P  _  jy  P 

-or  p  r  -a  r  ~  f 


and 


~~ag  K  +  i  P  in  P  Ay  _P 

-»t  r  f  t  ~  ?  -at' 


(ii. in 


II. 4.  NON-DIMEHSIONAL  VARIABLES  AND  PARAMETERS 

The  next  step  is  to  cast  the  governing  equations,  Eqs.  (II. 1)  -  (II. 3),  in 
their  most  convenient  non-dimensional  form.  The  various  non-dimensional  va¬ 
riables  and  parameters  can  be  properly  specified  by  referring  to  Fig. 11.1 

which  represents  the  various  blast  wave  coordinates  in  the  time-space  domain. 

« 

The  origin  of  the  system  is,  as  a  rule,  at  i:  -  O  and  IT  =  o  ,  and  the  co¬ 
ordinates  of  the  front  are  tn  and  .  With  reference  to  Fig.  n.i,  the  Eule- 
rian  space  profiles  are  obtained  by  seeking  a  solution  along  t  =  or 
C  =  i  ,  where 


(II.  i2) 


Therefore,  the  following  groups  of  non-dimensional  variables  and  parameters 


will  be  used: 
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i)  Physical  space  coordinates,  the  independent  variables  of  the  problem: 


(11.13) 


where  the  first  is  referred  to  as  the  field  coordinate  and  the  second  as 
the  front  coordinate.  The  symbol  (  \T  )  denotes  the  space  coordinate  of  a 
point  in  the  flow  field  while  (  ^  )  is  the  radius  of  the  front  ul  ihe 
same  instant  of  time,  subscript  (  O  )  specifying  a  reference  value.  It 
should  be  noted  that  when  the  front  trajectory  is  Known,  J  becomes  a  mea¬ 
sure  of  time. 

I 

iij  Front  parameters,  i.e.  variables  pertaining  uniquely  to  the  front 
motion  and,  therefore,  are  functions  of  only  the  frorit  cooroinate.  J  : 


d m  ^  ^ 

K1  " 


2.  ^  -  Jh  ii  _  2. 

TGj 


•  p  d (a 4*  _  J 
J  a 


(11.14) 


where  W„  *  is  the  front  propagation  velocity,  A  la  the  decay  co- 

m  %• 

efficient  which  expresses  the  deceleration  of  the  wave  front,  d  is  the 
velocity  o-p  sound,  f4  is  the  Mach  number  of  the  wave  front,  is  the  front 
velocity  modulus,  and  f  is  the  ambient  density  parameter,  while  subscript 
[2  )  denotes  conditions  of  the  ambient  atmosphere  into  which  (he  front  of 
the  blast  wav9  propagates. 

iii)  Gasdynamic  parameters  of  the  flow  field,  tse  dependent  variables 
describing  the  structure  of  the  flow  field: 


f  *  J£.  j  h  m  JL  J  9  v  a* d 

r  4  4  < 


(11.15) 
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iv)  Reduced  variables;  coordinates  of  the  phase  plane: 


Fm£  =  (£..«)  J  z  *  £3  =*,£.■) 

%  K<r  /l  >  J  V  ‘g) 


(11.16) 


v)  The  source  term  parameters: 


j  (a)  m  JL  -fi-  find  cJ  A.  JS_ 
n  fi  vn  F  r  **  £  ws 


(11.17) 


II. 5.  TRANSFORMATION  TO  NON-DIMENSIONAL  FORM 

1 .  Transformation  of  variables 

The  independent  variables  of  the  problem  in  non-dimensional  form  will  be 
X  *  %  (  r.  t  )  and  J  ■  J  (  tT ,  t  ) .  This,  consequently,  means  that 
r  *  IT  (  X  »  J  )  and  t  =  "t  (  X  ,  J  ),  and  with  the  aid  of  the  definite 
of  Eq.  (11.13),  one  has  for  the  non-dimensional  transformation 


-Q  -  -o#  -o  ~~a  J  -a 
*7  ^  ^ 


But  for  the  Eulerian  space  profiles,  •— —  =  O  .  Thus: 

o  \r 


~o  -  _ L.  ~o 


(11.10) 


~~o  %  _  ~~~a  J  3 _ 

~*t  ^  t  ^  J 

"a  #  ,  d  *0  .  +  Z?  -I  .  -dJS.  .  • 

-a  yn  dt  dt 


and 
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where 

_  1  /y 

~  -*r  ■  -*r„  r„ 

3lI  =  A.  =  A.  ,  7 
•'•n';  'n  J 

Thus 


3 _ 

"at 


Wrt 

»*' 


-q 


+ 


1 


(11.19) 


2.  Conservation  equations  in  non-dimensional  form 


Substitution  from  Eqs,  (11.13)  -  (11.19)  into  Eqs.  (II. 1)  -  (II. 3)  yield  the 
following  non-linear  partial  differential  equations  for  the  conservation  of 
mass,  momentum  and  energy,  respectively,  in  non-dimensional  form,  as  given 
by  Oppenheim  et  al.  (1971): 


wher- 
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(11.23) 


and 


+(JL  -  i )  2FZ _ £.  [/+  ^Jkk  +(JL-i)  5l|l L 1 

*Kx  ->tA%  k  L*  J 


~  ^  ~  f  ~  (s'  +  4 —  4r  ^ 


u> 


n 


(11.24) 


Equations  (II. 10)  and  (II. 11)  are  also,  respectively,  reduced  to 


-txt'  „  <&  +  1  ^  $  ~*]gJl 

p  A  A 


(IC.25) 


and 


Using  Eqs .  (11.25)  and  (11.265  to  eliminate  s'  f  om  Eq.  (11.24),  one  gets 


(11.27) 


where 


4>c  s£r|~  [*£■  -ftJr  - )  Hi] 


(11.28) 
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Equations  (11.20),  (11.22)  and  (11.27)  together  with  the  definitions  of  Eqs. 
(II. 21),  (11.23)  and  (11.28)  represent  the  governing  equations  in  their  ge¬ 
neral  form  without  restrictions  imposed  by  specific  equation  of  state  and 
including  any  possible  scuices  of  mass,  momentum  or  energy  that,  may  affect 
the  flow  field. 

3.  Reduced  autonomous  form 

The  governing  equations,  Eqs.  (11.20),  (11.22)  and  (11.27)  may  be  expressed 
in  terms  of  the  reduced  variables  defined  in  Eq.  (11.16).  When  these  varia¬ 
bles  are  substituted  in  Eqs.  (11.20),  (11.22)  and  (11.27)  yield,  respec¬ 
tively: 


(11.29) 


-JL 

X 


_  d5 

rr.  F 


(11.30) 


and 


4-  a  -r  [tM  $L  ^ 


(11.31) 

r 


The  above  equations  can  be  considerea  as  a  set  of  three  algebraic  equations 
for  three  unknowns  that  represent  the  logarithmic  gradients  with  respect  to 
X  of  f  ,  fj  and  $  ,  respectively.  Solving  them  algebraically  for  these 
gradients,  one  obtains  the  following  autonomous  form  of  the  blast  wave  equa¬ 
tions  : 
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As  a  consequence  of  the  definition  of  2  ,  one  has  also  from  the  above 
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where 
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(II. 48) 


II. 6  BOUNDARY  CONDITIONS 

In  addition  to  the  conservation  equations  and  the  equation  of  state,  the  for¬ 
mulation  of  a  blast  wave  problem  is  completed  by  the  boundary  conditions  which 
can  be  found  by  solving  simultaneously  the  equations  of  continuity,  momentum 
and  energy  across  the  shock  front.  Accordingly,  the  boundary  conditions  de¬ 
pend  on  the  problem  under  study.  Generally,  there  are  two  types  of  boundary 
conditions: 

i)  for  adiabatic  point  explosion  problem,  or  when  the  heat  of  reaction,  in 
the  case  of  reactive  medium,  is  distributed  inside  the  flow  field,  the  bounda¬ 
ry  conditions  or  the  Rankine-Hugoniot  relations  will  be  determined  as  follows: 

Consider  the  shown  control  volume,  Fig.  II. 2,  which  contains  discontinuity  in 
flow  variables  at  the  shock  front.  Applying  the  conservation  principles  of 
mass,  momentum  and  energy  between  states  1  and  2,  ahead  and  behind  the  wave 
front  immediately,  one  gets,  respectively 

%  =  4  ^ 
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and 
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Changing  to  a  frame  of  reference  in  which  the  shocK  wave  is  at  rest,  the 
above  equations  yield,  respectively 


(II. 49] 
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where  is  the  front  propagation  velocity,  and  subscripts  a  and  n  denote 

conditions  of  the  ambient  atmosphere  into  which  the  front  of  the  blast  wave 
propagates  and  conditions  immediately  behind  the  wave  front,  respectively. 

ii)  in  the  case  of  reactive  medium,  the  heat  released  may  be  applied  to 
the  blast  wave  at  its  front  due  to  laser  irradiation  or  chemical  reaction. 

The  boundary  conditions  in  this  case,  or  the  Hugoniot  relations,  are  determi¬ 
ned  exactly  as  above,  except  for  the  energy  equation  which  takes  the  form 

^  +  %  =  <?„  +  JL  +  Q&fMS?  (11.52) 


where  %  is  the  energy  deposited  at  the  front  per  unit  mass  of  the  medium. 
Defining  the  non-dimensional  chemical  energy  parameter 


(11.53) 
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Therefore,  the  boundary  conditions,  Eqs.  (11.49)  -  (11.52),  will  take  the 
following  non-dimensional  form,  respectively 


(11.54) 
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and 

*  ^a.  +  yr  h  +  jr  &  +  ^  ^ 
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(11.57) 


For  adiabatic  point  explosions  or  for  reactive  medium  with  distributed  heat 
release  inside  the  flow  field,  Eq.  (11.56)  1b  applied,  while  for  reactive  me¬ 
dium  with  instantaneous  heat  release  at  the  front,  Eq.  (11.57)  is  applied. 


II.  7  INTEGRAL  RELATIONS 

l 

The  integral  relations  for  a  blast  wave  are  the  expressions  of  the  principles 
of  global  mass,  momentum,  and  energy  conservation. 

I .  Mass  integral  aquation 

The  mass  integral  states  that,  at  any  moment,  the  amount  of  mass  distributed 
inside  the  flow  field  should  equal  that  engulfed  by  the  front  plus  that  of 
any  other  source  that  may  affect  the  flow  field.  Thus: 
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where  tf  j  is  the  total  mads  content  inside  the  flow  field,  while  Hj  is  a 
geometrical  factor,  defines)  by 

rij  s  'i  rJ  +  £  (L  l)  (J-2) 


is  th8  mass  of  the  ambient  medium  engulfed  by  the  front,  thus* 
.  *  j  ,  0  k* i  *n  hi 
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°  oJ  j+l 


with  the  aid  of  So.  (II. 14) 


O  p  , 
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thus 


(11.59) 


and 

H^is  the  mass  that  may  be  added  to  the  flrw  field  by  a  source  from  the 
instant  of  explosion  till  the  instant  under  consideration,  thus 


ss  j  °  j  ff*  VI ^  Jr  Jr 

Q  4 


(11.60) 


Therefore,  the  substitution  of  Eqs.  (11.03)  and  (11.60)  into  Eq.  (11.06) 
yields  the  mass  integral 
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which  for  a  sourceless  flow  and  uniform  ambient  density,  =  constant,  re¬ 
duces  to  the  simple  form 


(Ui) 


(II. 62) 


On  the  other  hand,  the  mass  integral  may  be  obtained  directly  by  integrating 
Eq.  (II.  1),  the  mass  conservation  equation,  with  respec*  to  \T  ,  from  the 
center  of  symmetry  to  the  shock  front,  to  yield 


(fr1) 


r0‘  =  Jr 


(II .03) 


By  utilizing  Leibnitz  rule,  the  first  term  gives 

r*  *  ■ 

<*r  =  jtj''{‘rUr-?„  r„* 

which  may  be  substitute  1  into  Eq.  (11.63)  to  yield,  noting  that  dji.  «  . 

dt 

y*  »  K 

f  A  r  ■#*  fn  *n  ^  )  s  J 

But,  from  Eq.  (11.49),  one  has 
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The  last  equation  may  now  ba  integrated  with  respect  to  time  tu  ytOid 
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which,  after  integrating  the  first  term  on  the  right  hand  siae  by  part,  re¬ 
duces  to  the  required  mass  integral,  exactly  as  given  by  Eq.  til. 61). 

i. 

Now.  for  fa  =  constant,  and  noting  the  definitions  of  Eq.  (II. 14), and 
t dff  ,  Eq.  (11.17),  Eq.  (11.61)  may  be  non-dimensionalized  to  yield 

i 
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which,  in  turn,  for  a  sourceless  flow  reduces  to 

,  f'kt'Jx  « -4_ 
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(11.65) 


2 .  Nomen'. um  integral  equation 

The  momentum  integral  may  be  obtained  by  integrating  Eq.  (II. 2),  th«  momentum 
conservation  equation,  with  respect  to  f"  ,  from  f  -  o  to  K  1  ^  ,  as 

follows 
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The  first  term  yields.  after  utilizing  Leibnitz  rule. 
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Substituting  the  above  equation  in  Eq.  (11.66),  -o  ield 
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By  solving  Eqs.  (11.49)  and  (11.50)  simultaneously,  one  obtains 
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Integrating  the  above  equation  with  respect  to  time,  one  gets 
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(11.67) 


where  is  the  momentum  of  blast  wave. 

'i 

Hence,  for  constant  f  Eq.  (11.67)  may  be  non-dimensionalizt;d  to  give 
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3.  Energy  integral  equation: 

The  energy  integral  states  tnat  the  total  energy  of  a  blast  wave,  »  should 
equal  the  summation  of  the  energy  deposited  initially  to  generate  the  wave, 

£.  ,  the  ambient  internal  energy  engulfed  by  the  Front,  £  ,  the  chemical 

energy  liberated  in  the  case  of  reactive  medium,  ,  and  the  energy  due  to 

sources  within  the  wave,  /T  .  Thus,  for  JL  =  constant: 
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where 
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if  the  heat  liberated  is  considered  at  the  front,  thus 


U  i 


while,  in  the  case  of  distributed  heat  release,  £ck 
where 


will  be  obtained  as  £  , 


Therefore  the  energy  integral  will  be,  in  general 
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For  a  sourceless  flow  and  heat  released  at  the  front,  the  energy  integral 
reduces  to 


ni 


,  i&vhi 
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(11.71) 


We  are  going  to  use  the  definition  of  Y"0  ,  the  reference  radius,  as  that  gi¬ 
ven  by  Oppenheim  et  al.  (1971)  which  includes  the  two  most  important  initial 
parameters  (  £;  and  )  of  the  problem: 
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Using  the  deiinition  of  to  replace  Ej  in  Eq.  (11.70)  by 


£ 


ro 


Thus,  the  energy  integral  equation  in  its  general  form,  Eq.  (11.70),  can  be 
expressed  in  non-dimensional  form  as  follows 


which,  for  a  sourceless  flow  field,  reduces  to 
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where 


<5~  = 

*  W/ 


=  id  J&. 


Differentiating  the  above  equation  with  respect  to  *d  and  using  the  definition 
of  A  *  one  can  easily  obtain  an  expression  for  the  decay  coefficient 


A  _  0+*)  ^3  -  -  %  d 

-  ij  dp's 

jy 


(11,75) 


wtiile,  for  adiabatic  point  explosion,  the  above  equation  reduces  to 
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11.8  CASE  OF  PERFECT  GASES 

As  formulated  above,  the  introduction  of  a  general  equation  of  stat e,  e(P,f)» 

leads  to  two  thermodynamic  properties,  namely  the  internal  energy  and  the 
speed  of  sound  factor  ,  which  include  essentially  the  characteristics  of 
real  gases.  However,  due  to  the  complexity  of  the  governing  equations  in 
their  general  form,  the  assumption  that  the  flow  field  behaves  as  a  perfect 
gas  with  mean  specific  heat  ratio  has  a  wide  popularity  in  most  literature. 

i 

1 .  Properties  of  perfect  gases 
The  equation  of  state  of  a  perfect  gas  is 

P  =  f  RT  (11.77) 

where  R,  is  the  gas  constant  and  f  is  the  temperature.  One  has  also  the  fol- 
liwing  relations: 
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3.  Boundary  conditions 

Equation  (11.54)  of  the  boundary  conditions  remains  unchanged,  while  Eqs. 

i 

(11.55)  -  (11.57)  reduce,  respectively,  to 


(11.85) 
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For  any  given  shock  strength  d  ,  there  exist  two  possible  solutions  repre¬ 
sented  by  the  two  signs  in  Eq.  (11.91).  The  plus  sign  denotes  the  overdri¬ 
ven  detonation  solution,  and  the  minus  sign  the  weak  detonation  solution. 

The  unique  Chapman-Jouguet  solution  is  obtained  when  the  two  solutions  coin¬ 
cide,  i.e.  when  the  terms  between  square  brackets  in  Eq.  (11.91)  equal  zero, 
yielding 


y  r  =  [<*t  o  %  *  0  ±  /[tfriD.f  -  ifZT 


(11.92) 


with  tne  minus  sign  applying  for  detonation,  while  the  plus  one  for  deflagra¬ 
tion.  In  general,  only  the  overdriven  detonation  solutions  are  of  interest; 
and  one  may  discard  the  bottom  sign  in  Eq.  (11.91)  for  weak  detonations. 

Thus; 
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The  other  two  boundary  conditions  for  ^  and  are  still  applied  through 
Eqs .  (11.54)  and  (11.65). 

4.  Integral  relations 

The  mass  and  momentum  integrals  in  the  case  of  perfect  gases  are  not  changed 

and  have  the  same  forms,  given  by  Eqs.  (II  64)  and  (11.68),  except  for  the 

momentum  integral,  Eq.  (11.68),  P  changes  to  %  .  While  the  energy  inte¬ 
nt 

gral,  Eq.  (11.73),  reduces  to 
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which,  for  a  sourceless  flow  field,  becomes 
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and,  the  decay  coefficient  takes  the  form 
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which,  for  adiabatic  point  explosions,  reduces  to 
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II. g  CASE  OF  SELF-SIMILARITY 

The  ^e If -similar  motion  of  a  medium  is  one  in  which  the  number  of  independent 
variables  in  the  fundamental  system  of  equations  is  reduced  from  two  to  one, 
thus,  the  system  of  non-linear  partial  differential  equations  is  reduced  to  a 
system  of  ordinary  differential  equations.  However,  the  value  of  this  solution 
is  confined  to  the  early  time  regime  when  the  shock  wave  is  strong  enough  to 
neglect  the  effect  of  counter-pressure.  This,  consequently,  means  that  the 


-  30  - 


front  propagation  velocity  is  very  large  approaching  infinity.  Thus, 
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(11.99) 


1  *  Governing  equations 

Sutistituting  this  condition  of  self-similarity  into  governing  equations,  Eqs. 
(11.20),  (11.22)  and  (11.27)  yield  respectively 
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(11.102) 


For  a  perfect  gas,  the  energy  equation  reduces  to 


(11.103) 


where  is  given  by  Eq.  (11.84). 

Thu  above  set  of  equations  is  an  ordinary  differential  one,  which  can  be  solved 
numerically  once  the  source  terms  are  specified. 

j 

I 

2 .  Boundary  conditions 

The  boundary  conditions  in  this  case  are  simply  obtained  by  putting  id  =  O  in 
Eqs.  (11.55)  -  (11.57),  while  Eq.  (11.54)  is  still  applied  without  any  restric- 


tions.  Thus,  one  has  the  following  relations 


(II. 104) 


(11.105) 


where  Eqs.  (11.56)  and  (11.57)  reduce  to  a  single  equation  given  by  Eq.  (11.105) 
This  means  that  the  heat  released  at  the  front  is  negligible  in  the  self-simi¬ 
lar  solution  compared  to  the  initiation  energy,  ,  which  is  the  predominant 
parameter  governing  the  flew  field. 

For  a  perfect  gas,  the  boundary  conditions,  Eqs.  (11.54),  (11.104)  and  (11.105), 
yield  the  following  relations: 
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3.  Integral  relations 

The  mass  integral,  in  this  case,  is  given  by  Eq.  (11.65),  since  for  strong  ex= 
plosions  the  double  integration  in  the  right  side  of  Eq.  (11.64),  which  includes 
the  effect  of  mass  source  term,  will  equal  to  zero.  On  the  other  hand,  the 
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energy  integral  still  takes  the  form 
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(11.109) 


while  the  right  side  or  Eq.  (11.73)  leads  to  a  trivial  solution  of  cTj  =  O  . 
Thus,  the  energy  integral  is  obtained  by  numerical  integration  of  Eq.  (11.109) 
from  the  flow  field  solution. 

liquation  (11.109)  is  reduced  to  the  following  relation  in  the  case  of  perfect 
gases 
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while  the  decay  coefficient  takes  the  following  simple  form 

X  ~  4+1 


(II. Ill) 


H.10.  COMPATIBILITY  EQUATION 

At  any  time,  the  particle  velocity  at  the  center  of  symmetry  must  be  equal  to 
zero,  that  is 


u  (o ,  t)=  o 

This  represents  a  compatibility  equation  which  must  be  satisfied  by  the  correct 
solution.  This  equation  in  non-dimensional  form  can  be  written  as: 


f  (OjH)  =  o 


(11.112) 


Figure  Captions 


Fig.  II. 1.  Blast  wave  coordinates 

Fig.  II. 2.  Control  volume  of  a  blast  wave 


CHAPTER  III 


THE  QUASI-SIMILAR  SOLUTION 


III  .1 .  INTRODUCTION 

Thi  conservation  equations,  as  derived  previously,  are  highly  non-linear,  r.on- 
homogencous,  coupled,  partial  differential  equations.  Therefore,  a  closed  form 
solution  is  so  far  impossible,  and  approximate  nolutions  appear  tn  be  essential. 
Mainly,  three  analytical  techniques  are  used  for  predicting,  the  non-self-similar 
flow  field  of  a  blast  wave,  the  perturbation  method  of  Sakurai  (1965),  the  qua¬ 
si  similar  method  of  Oshima  (1960,  1964)  and  the  power  density  law  of  Bach  and 
Leo  (1970). 

In  this  chapter,  the  quasi-similar  technique  will  be  discussed  briefly  with  a 
complete  description  and  derivation  of  this  method  of  solution  to  be  applicable 
for  any  type  of  problems. 


III. 2.  TAYLOR'S  SERIES  EXPANSION 

To  obtain  the  approximate  quasi-similar  solution,  one  has  to  perform,  firstly, 
a  double  Taylor's  series  expansion  on  the  dependent  variables.  As  mentioned 
previously,  any  dependent  variable  of  the  problem  is  a  function  of  the  two  inde¬ 
pendent  variables  X  and  i)  .  Thus,  for  any  dependent  variable  5^  ,  the  appli¬ 
cation  of  the  double  Ta>lur's  series  expansion  yields 
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where  o.  and  b  are  the  values  cf  X  and  y  ,  respectively,  in  the  neighborhood 
of  which  the  solution  is  to  be  obtained,  provided  that  the  function  (Xj*!) 
together  with  its  partial  derivatives  up  to  1  lie  order  are  continuous  in  the 
neighborhood  of  the  point  (  a  ,  b  ). 

Since  one's  interest  centers  upon  the  region  of  the  flow  field  near  the  wave 
front,  for  any  value  of  Mach  number,  onu  may  set 


a.  =  i 


b  =  y 


Eq.  (II1.1)  then  reduces  to 
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It  now  becomes  possible  to  substitute  the  above  equation,  of  course  with  7  -  f  , 
h  .  and  9  ,  in  the  equations  of  motion  and  proceed  to  obtain  the  expansion, 

soJution  in  the  normal  fashion. 

With  reference  to  the  autonomous  equations  (IJ.32)  -  (11.34),  only  the  quanti¬ 
ties  subscripted  by  J  will  be  expanded,  thus  yielding  a  set  of  first-order 
ordinary  differential  equations  identical  to  those  of  the  self -similar  case  ex¬ 
cept  that  the  dependence  of  J  ,  or  y  ,  is  included  only  as  additive  algebraic 
terms.  For  any  dependent  variable  ,  the  autonomous  equations  give 


o 


fill. 3) 


where  0  *  and  &  are  defined  in  Eqs.  (11.36)  -  (11.48) 

pc  J  tv 

equation  may  now  be  put  in  the  form 


In  <f  _  <7#  +  -t  $1} 

***  xu 


wl 


iich,  by  differentiating  with  respect  to  A  J  ,  yields 


-  ( &  \  s  ^ _ r  f*  +  &■}+ &  i 

*  K  ~*Pn  J  ~  L  ^  C>  J 


where  the  order  of  differentiation  of  the  term  on  the  right  hand  side 
interchanged.  Thus,  with  the  aid  of  Eq.  (III. 2),  one  has 


zAJL  *(3A£)  fa — I  +, 

->M  -»M  *.<  L-'A.j  a-t 


•f  —  — 


or 


=(iA£ 


CxA£\  -({.*>[ 2;-  (5ul3  ±Ll)  1  * - 


1  /n  A  4A 


M 


-C'A^1)  -  0 -%)  f /  Hr  -r  +  //qi  )1 


+ - - 


and 


llal  — .  f(gx  +  <*]  *  Gj?'1  l~l  t _ 

•4j  *o  Jlsi 


where,  denoting 


•U 


The  above 


(III. 4) 


has  been 


l III. 8) 


lill.5u) 


(III. 5b) 


(III.5c) 


by  a  prime, 


(III. 6a) 


V 


3 _  (  %  +  F3  +  &  ) 

^  %D 


and 


-ff*-/0?— >  «&><<♦«!♦ 


'  ^  =  f  VP-  ifi)  +  fF.J)  (F'_  |1  ) 

+  (rr'+Fr')-£] 

.J-(Z.«U)rF-£][z'-  Z(f-'.f)] 


f.i 


>  r 


t 

rF$ 


u 


3'F  z'-z(r',£+r:)i  1 


(III. 6b) 


(III. 6c) 


(III. 7) 


(III. 8) 
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fJ  rn  =(±-n  4>'f  -  ?  % 

~*Sn J  r  r  ,  r 


while 


(III. 95 


<  •  5ft  -  &>  K*  F*  F f«] 

*  D[(4ij  f'+  £"]  +  ^  -£l 


(i-r) 


(III. 10) 


(III. 11) 


N /  s  s  -J —  fFF^*  f'f  -  #L  -  &'  D1 

w  rf-nl  ^  •*>  tf  ■* 


a 1  <*-F> 


f/  yJEl) 

^  1  i-r 


mi. 12) 


<  -  5ft-  =  7^5  f  r'F*  *  r«'>  *  r>S  'lrF<U» 


-A*/]  D/+D[(j+i)(rF/i-r'F)-fi'+£]  +  &JiL  (III. 13) 

L  1  -*  r 
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/ 


1 _ 

<*-n 


{ nrtj  +  r/j' ) + rV  fj 


+ 


(III. 14) 


(III. 15) 

where 


2fl-F)F'] 


(III. 16) 


Now,  for  simplicity,  consider  the  case  of  perfect  gas  with  constant  specific 
heats  and  constant  f  .  For  adiabatic  point  explosion,  the  boundary  condi¬ 
tions  are  given  bv  Eqs.  (11.88)  -  (11.90),  and  since 

r'm  zJL  = 

one  may  then  evaluate  Eqs.  [III. 7)  -  [III. 161  at  X  =  1  ,  by  noting  that 


.  JS_  (i.H) 


lf*f  * 


) 
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<**><*  +  *#> 


h  =  - - 

n  v-  i  +ziJ 


la'  B  (M')  JziBtii L 

rt  fy-lai3)‘ 


h"  -  +  a  +  y  #] 

""  (y-(+4!0*  l-y.ia.i3  JyJ 


^  y-  (y-1)^ 
V  (*•+•  1) 


<  = 


II 

%* 


F-‘L(»)  (A 


Ai  ) 


z  -  -  .J-  [A  ir  -  car- 4)31  fir- l  +iy] 

n  “  k*  Tivi?  L  J  L  J 


Z*  =  -AM-  f  4  [£*-  ('8'-!)^]  -  (**-*)  («■-!  +  *3)1 

<*vl/  L  J 


The  above  equations,  then,  give,  for  any  order  of  the  expansion,  algebraic  ex¬ 
pressions  for  2^  .  .  and  .  These  may  be  substituted 

in  /v  ,  H j  ,  and  »  respectively,  to  yield  a  set  of  first-order  autono¬ 
mous  differential  equations  with  respect  to  the  field  coordinate  %  ,  where  ^  , 
%  and  its  derivatives  appear  as  parameters.  The  gasdynamic  flow  fields  may 
now  be  obtained  for  each  value  of  jj  .  where  ^  and  its  derivatives  are  found 
either  from  experiment  or  from  the  energy  integral,  by  integrating  these  auto- 
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nomous  equations. 


III. 3.  QUASI-SIMILAR  METHOD  OF  SOLUTION 

The  concept  of  the  so-called  Oshima’s  "quasi-similar”  solution  states  that  the 
derivatives  of  all  dependent  variables  with  respect  to  the  front  coordinate, 

J  or  Ij  ,  are  considered  constants  in  all  of  the  flow  field  and  equal  to 
those  corresponding  to  X  =  i  .  This  concept  arises  from  the  fact  that  the 
variations  of  the  dependent  variables  f  ,  ji  ,  and  9  are  not  sensitive  to 
the  independent  variable  J  ,  or  ij  ,  as  compared  to  their  variations  with  the 
other  independent  variable  %  .  However,  this  assumption  causes  the  solution 
to  be  accurate  for  high  Mach  numbers,  corresponding  to  small  values  of  9  , 
and  just  behind  the  wave  front. 

From  the  mathematical  point  of  view,  this  method  of  solution  may  be  considered 
as  the  zeroth  order  of  the  Taylor’s  series  expansion  of  a  function  of  two  in¬ 
dependent  variables.  It  follows  from  Eqs.  [III. 5)  that 


iL 


-  p£) 


(III. 17a) 


x=i 


(III. 17b) 


*>9  \ 

xGi 


(III.17c) 


This  leads  Oshima  to  express  the  non-dimensional  dependent  variables  in  the 
following  form 
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=  %<*)  ■£(*) 

k(Xj2)  =  <*) .  k^) 

S(X,d)=  ^(x'i.^Cd) 


Differentiating  the  above  equations  with  respect  to  ,  one  obtains 


m  p(X*H)  d  £ 

^  y  £  a)  '  d  y 


skcXjj)  =  k(X^)  J  kj  CJ) 

“  k^ft)  '  <1* 


-a  9(aotJ)  _  <3(x,  *)  .  d  9*  C*) 
^  "*  %  CV  '  d  a 


(III.18a) 


(111.18b) 


(111.18c) 


Taking  the  values  of  Px  (  9  ),  ki  (  9  1  and  9*  (  9)  as  their  values  at  the 
wave  front,  yields  to  the  required  solution. 

For  constant  and  when  the  medium  may  be  considered  as  a  perfect  gas  with 
constant  specific  heats,  Eqs.  (III. 18)  may  be  obtained  simply  from  the  boundary 
conditions  of  adiabatic  point  explosion,  Eqs.  (11.88)  -  (11.90),  or  from  the 
boundary  conditions  of  reactive  detonating  blast  waves,  Eqs.  (11.54),  (11.85) 
unil  (11.93). 

1 .  Adiabatic  point  explosion 
In  this  case,  one  obtains  from  Eqs.  (11.88)  -  (11.90) 


df.QQ 

dtl 


-2 

&«■! 


( III . 19a) 


d  * 


A  A»  c^d)  _  d  htfi 

"  djf 


a  (y-t-i) 
(*-i  +*9)a 


(III. 19b) 


d  %  «Q  __  J  9p  _  (y-  1) 

d9  ~ 


(III.19c) 


Substituting  Eqs.  (11.60)  -  (11.90)  as  well  as  Eqs.  (III. 19)  into  Eqs.  (III. 18) 
yields,  respectively 


2L  -  P 

(i-tf) 


(III. 20) 


~*k  s  _  2  h 

and 

-a  9  _  (x-i)jl _ 

~>d  zv-v-w 


(III. 21) 


(III. 22) 


Substitution  of  Eqs.  (III. 20)  -  (III. 22)  into  governing  equations,  Eqs.  (11.20), 
UI.22)  and  (11.83),  yields,  respectively 


jfzjt)  Af  i  A9 

f  Ax  kf  Ax 


xe 


(III. 23) 


(111,24) 
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(P-X)  J9 

3  *iX 


_  (JL-Xl  Ah.  +  Ac  8  <f> 

K.  Ax  E 


where 


(III. 25) 


/?  H 

(V-l)  +  X  1 


0  s  - . .  y  ±-i— 

*  (i-tO 


G  3  -  [ _ ll - : -  4.  Xfi  ] 


Thus,  the  quasi-similarity  assumption  has  reduced  the  set  of  non-linear  par¬ 
tial  differentail  equations  to  a  set  of  ordinary  differential  equations  with¬ 
out  making  the  strong  shock  approximation. 

Solving  for  the  derivatives  Ah  ,  Ah.  and  •—  in  Eqs.  (III. 23)  -  (III. 25), 

Ax  AX 

one  obtains 


J£  S  igg  +  ^g£  +  k<#r  -  %d)  -  V3(#rt  -I/O  -  s 

k(P-x)*-X3 


Ml  *  f  AL  +JJL  +  A/?  -  *&J 

Ax  (f-x) 


=,  ./t  [<P-x)j£  +  AZ6  -Z*^] 


(III. 26) 


(III. 27) 


(III. 26) 
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2.  Detonating  blast  waves 

Applying  the  same  previous  procedure,  one  arrives  finally  at  the  same  Eqs. 
(III. 26)  -  (III. 28),  except  for  the  parameters  f\  ,  Q  and  C  ,  which  will, 
with  the  aid  of  Eqs.  (11.54),  (11.85)  and  (11.93),  take  the  following  form 

/?  »  .  g — 
i-f 

rn 


E>  s.  SB-  —  q.5 

& 


a  e  *<r±3..  - 1 .  y* 
*1 Z  +  * 


where 


m.4.  ADIABATIC  INTEGRAL 

In  order  to  find  a  check  on  the  numerical  method  used  for  solving  the  governing 
equations  of  the  quasi-similar  solutions  in  the  case  of  sourceless  adiabatic 
flow  fields,  one  may  derive  an  intermediate  integral,  the  so-called  adiabatic 
integral.  The  continuity  and  energy  equations,  Eqs.  (HI. 23)  and  (III. 25), 
can  bo  rewritten  for  sourceless  flow  fields  as 


<f-Z) 


Jjnh  Jf 


i  X.  +  A  A  x  o 


(III. 29) 
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(P~oc)  dJkl  -  y  (f-x)  dJuk  +  /Xc  Go  chi. 30) 

Ax  Ax 

with 

A  -  -  -  11 _ 

^—1  +  x  $ 


and 


c  «  -4*  r  — 1 - 1 _ i 

L  ay  -«r-i)y  y-l  J 


To  fipd  the  adiabatic  integral  from  the  above  two  equations,  one  may  multiply 
the  first  of  them  by  an  integration  factor  y'  and  then  subtract  Eq.  (III. 30) 
from  it,  one  obtains 

(f-x)  O'*  m  _  A?ij±]  +(d£  .  i)  f  ♦  >£  v 

+  &  +  (V'A-  c  )  mo  (HI. 31) 


For  the  above  equation  to  be  integrated,  the  value  of  the  integration  factor 
¥  will  be  chosen  as  to  satisfy  the  condition 

JV  =  -  \  C*  A  -  C  )  -  * 

< 

from  which  one  gets 

Sft  _  X 

1+  4  +  CX  A 


Or,  with  the  use  cf  the  definitions  of  >  and  G  , 
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✓  *  2  [tf/or-i+atf) -i /(**»*- ya)3 

i+i  -  /  (*.  i  +*ti) 


(III. 32) 


Therefore,  Eq.  ( III . 31 J  will  be 


(¥'*.  *)  iiil  >  l{o3.  +  '*'($LzJ-)  +  Jjt! 

dX  Kf.X  * 


s  O 


which  yields  upon  integration 


[A  xJ(x-P)]  s  K 


where  K  is  a  constant  of  integration  which  can  be  determined  from  the  boundary 
conditions,  cqs.  (11.88)  -  (11.90),  one  gets 


J® 


K  .  ii_  <l-£) 

“(1 


y  <r-»i)  [ur+i)/<r-i  +*■*)[ 

(a  if-  Xd  +  d) 


Therefore,  the  adiabatic  integral  will  be 


zV\l*/Qr-  i+  2i)  -  * 

Jhf  x\x -f)  J  i+i  ~2.'k'A/(s-i  +  i'4)  -  (III. 33) 


,  d  =  0  and  X  =  -*+  i  ,  the  above  integral  re¬ 


in  the  case  of  strong  explosions 
duces  to 


1 
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f- 

I 

I  J 

I- 

>  , 


» 

s 


a!  [kx'ix-f)} 


=  (*$*-)  <■ 


(III. 34) 


The  adiabatic  integral,  Eq.  (III. 33),  has  been  used  as  a  check  on  the  precision 
of  machine  computations.  The  step  size  A  X  in  the  Runge-Kutta  method  was 
assigned  different  values  until  the  values  of  f  ,  k  and  $  ,  obtained  from 
the  numerical  calculations,  satisfied  Eq.  (III. 33)  in  a  desired  accuracy. 


* 

! 

i 


m.5.  SINGULARITY  ANALYSIS 

In  order  to  analyse  the  singularities  of  the  governing  equations,  it  if  more 
convenient  to  transform  them  into  the  pihase  plane,  i.e.  in  terms  of  the  reduced 
variables  f  and  2  defined  by  Eq.  (III. 16),  which  are,  for  a  perfect  gas 

2(40*)  =  and  RXj*)  =  (III. 35) 

X'h  K 


In  the  governing  equations  for  quasi-similar  solutions,  one  has  derivatives 
with  respect  to  fc  only  while  the  variable  y  appears  as  a  parameter.  There¬ 
fore,  one  can  write 

-i_  4*  =  i  Al  _i_A  .  x. 

z  Jx  9  3*  Id*  X 


and 


JLi£ 

r  mc 


i  t 

T  * 


The  above  two  equations  combined  with  Eq.  (III. 35)  when  substituted  in  Eqs. 
(III. 23)  through  (III. 25)  give  for  a  sourceless  flow  field 
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<n  i )  AMl  +  r  JA *  (u  i  +  JA) 

dr  Af  r 


+  i 


-f-  i&A  +  jja£  [F(F-  i)  +  JLZ.  +  IFB]  +JL  +(F-i)mO 


(i-K)(F„i)  J&A  +  A!s£.  (*(r.i)  +  Xc\  +  (F-i)  Afi*  ,.0 

dr  af  «■  j  J r 

Solving  for  the  logarithmic  derivatives  A /*!.A  ,  d anCj  AJ?!L 

Ar  Af  af 

gets 

Af .A  .  if-dUai+Jf)  -,\F(F-i)&  +  Z(XC/t) 
df  (r.i)  {  F(F- *)*+  X  F(F-t)  0  -  Z[rti+i)  +  A%i-  W*] 


A  fox  - _ £  -  ^ _ 

dr  F(F-i)i+XF(F-i)Q-7  [f  U+i)  +  A%  +  XC/rf 


A(oZ  ,  -A  F(F-i)Q  »  ZAcT/d].^y.j)^c3^-ft.f?‘ 

dF  <se-i){F(F-±)\%F<F  0  8-  Z  [rfiVi)  +  /4A  +  %£]} 


(III. 36} 


(III. 37) 


(III. 36) 


,  one 


(III. 39) 


(III. 40) 


(III. 41) 


The  Rankine-Hugoniot  conditions  at  the  shock  front  can  be  rewritten  in  terms 


of  the  variables  Z  and  f  by  combining  Eqs.  (11.80)  through  (11.90)  with  Eq. 
(11.16) 


fa*  -  Qf--i)'d]  (*-*  +  *&) 


(III. 42) 


-  Fa**)  =(S-l)  <*- U)  (III. 43) 

o +  * 

Eliminating  H  from  the  above  two  aquations,  one  obtains  a  relation  between 
Zn  and  .  It  is 

Zr\  =  -  Z\  )  (*  +  fn  )  (III. 44) 


It  is  obvious  that  the  solution  of  Eqs.  (III. 26)  through  (III. 28)  will  diverge 
when  the  denominators  vanish.  This  will  happen  at  the  values  of  X  given  by: 

i 

X-O  (III. 45a) 


x  *F 


(III. 45b) 


2  * 


(III.45c) 


The  divergence  of  the  solution  at  X  =  0  may  be  avoided  by  expanding  the  fluid 
properties  in  terms  of  X  near  the  center  of  symmetry.  The  continuity  ariu 
momentum  equations  for  a  sourceless  flow  field,  after  some  algebraic  manipula- 
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tions,  can  be  written  as 


_ 


=  o 


(III. 46) 


(Pa1)  +  ^  =o  (III. 47) 

Integrating  Eq.  (III. 46)  with  respect  to  |r  ,  one  obtains 

*  \r  *' 

fur*  s  -  /  r*  Ay-  ♦  Cot ut.. 

oJ  -at  * 

undef  the  condition  that  U  a  0  at  r  *  0,  the  constant  in  the  above  equation 
will  be  zero,  and  this  equation  becomes 

U(r''*)  s  -j?rJ r*  (?f>  d'r  (in. 46) 

In  a  similar  manner,  Eq.  (III. 47)  upon  integration  with  respect  to  r  ,  yields 

p(r,t)  =  -  fUz  -  f  Av~-f±.  fa*  At-  +  G>*rt. 

0J  f  r  J  a 

Also,  in  order  to  determine  the  constant  in  the  above  equation,  we  use  the 
* 

condition  that  U  =  0  at  y  -  0.  Hence,  the  above  equation  becomes 

P<r,t>  =  P(o,t)  -  fu3-.  /  ?Lfui  At-  -  fim1  At- 

«  "o  t  9J  r 


(III. 49) 
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Differentiating  Eq.  (III. 48)  with  respect  to  time  and  substituting  into  the 
above  equation,  one  gets 

P(r/t)  «r  P(o,t )  -  fu*-  ij  Jir  +/-4p  (HI. 50) 

*  •  *• 

Assuming  that  the  density  can  be  expanded  in  the  following  power 

series  near  the  center  of  symmetry 

P(iot)  s  r*  [  fQ  <t>  +  fx(U  r  +  f2  rl+ - ]  (III. 51) 

where  ot  is  a  constant  and  the  coefficients  £  ,  j •  fx  ,  ~  are  functions 
of  time  only. 

Substituting  the  proceeding  series  into  Eq.  (III. 48)  and  integrating  with 
respect  to  K*  ,  one  obtains  the  velocity  profile  represented  as  a  series  of 
the  form 


a<r,t>  e  r  +  Utct)  r  +  ux(t)  (III. 52) 

i 

where  the  coefficients  U9  ,  Ut  •  , --  are  functions  of  time  only. 

Substituting  from  Eqs.  (III. 51)  and  (III. 52)  into  Eq.  (III. 50),  one  finds 
that  the  pressure  can  also  be  represented  as  a  series  of  the  form 

=  Pco,t)  +  c-t)  +  p,( i)  Px(t)  rz.+  __J  (in. 53) 

where  P(ojt)  >  P0  *  Pt  •  %  _ are  functions  of  the  time  only. 

For  convenience,  the  three  series  for  the  density,  velocity  and  pressure  can 
be  rewritten  in  a  non-dimensional  form  as 


h.  m  pc  ( A*  +  A/  9c  +  A*  ^  ^ 


(III. 54) 


f  =  -f  ^  ^  *** - )  (III. 55) 

Je  <3(0,*)  +  (%  +  9±X  +  3tKl+—  )  (III. 56) 

where  the  coafficienta  A#  .  A*  .  A*  *  ---*  j£  •  4  *  4  »-—  *nd 
3(0^1/)  '  9®  *  9i  .  9a  • - ara  functions  of  only. 

It  is  evident  that  for  small  values  of  X  »  i.a.,  x  « 1 •  the  following 
asymototic  formulae  hold 

f  *  f  X  (III. 57) 

O 

kxk4X*  (HI.  56) 

9  r  9(o, a)  4.  9.  ****  (III. 50) 


Tha  fact  that  the  velocity  profile  is  linear  near  the  center  of  symmetry  in¬ 
dicates  that  the  value  of  PM  will  be  finite  and  no  divergence  of  the  solu¬ 
tion  will  occur. 

Because  one  cannot  proceed  with  the  numerical  solution  till  tha  cantor  of 
symmetry,  the  solution  will  be  based  on  matching  the  numerical  aulutiun, 
from  X  «  1  to  X  •  32  ,  with  the  asymptotic  formulae  fur  f  ,  A  and  9 
given  by  Eqs.  (III. 57)  -  (III. 59),  provided  that  the  value  of  %  must  be  as 
small  as  possible. 

In  terms  of  the  reduced  variables  f  and  Z  ,  the  other  two  singularities, 
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The  value  of  the  variable  f  at  the  shock  front  is  given  from  Eq.  (III. 43) 
while  its  value  at  the  center  of  symmetry  is  obtained  by  substituting  X  B  0 
in  Eq.  (III. 55).  It  gives 


-4) 

X  %mo 


Therefore,  the  range  of  variation  of  F  is  given  by  the  inequality 

l 

P.  4  F  « 

The  relationship  between  Z  and  F  has  been  obtained  by  solving  Eq.  (III. 41) 
numerically  as  F  ranged  over  its  value  at  X  -  1  to  its  value  at  X  ■  0. 

Since  the  value  of  ^  is  not  known  a  priori,  an  iteration  procedure  is  per¬ 
formed  to  search  for  the  correct  A  •  The  correct  %  may  be  obtained  based 
on  satisfaction  of  the  energy  integral  equation,  Eq.  (11.96),  or  it  may  be 
found  when  it  satisfies  the  condition  of  zero  velocity  at  the  center  of  symme¬ 
try,  Eq.  (11.112).  When  the  value  of  is  greater  or  less  than  the  correct 

one,  the  solution  will  diverge  at  some  non-zero  values  of  X  corresponding  to 
the  singularity  f  =  1  or  Z  s(i-F)X 

Furthermore,  from  the  theory  of  ordinary  differential  equations,  the  slope 

jf:  will  take  the  indetermined  value  when  both  the  numerator  and  denomi¬ 

nator  on  the  right  hand  side  of  Eq.  (III. 41)  vanish  simultaneously.  This  will 
correspond  to  a  number  of  singular  points  given  as  the  intersection  of  the 
Isocline  of  zero  slops  with  the  isocline  of  infinite  slope. 


The  isocline  of  zero  slope  is  given  from  the  following  equation 


[jr(r-i)  -AFs]  .zz  +  z  (F~i)x +  <r-4)Ac}]  =o  im.t 


while  the  equation  of  the  isocline  of  infinite  slope  is 


(F.i){r(F~i)  +  'XF(F-i)  8  -  zr<ui) -z  %(*  +  £/*)]  (in. 63) 


Solving  for  Z.  from  Eq.  (III. 63),  one  gets 

Z  =  F(F-i)%  +  ^F(F-i)3 

F  <j+i)  +'X(A+  C/H') 

Substituting  the  value  of  Z  from  the  above  equation,  one  obtains  an  expression 
for  the  values  of  F  at  the  singular  points.  These  values  are  obtained  from 
t  he  roots  of  the  equation 

(F-i)  [  or-  i)  {  iF(F :  i)  -  A  FB  }  +  2  (F.  if  +  A  (F-  1)  C  J 


F(F-i)x+  %F(F- 1)  8 

F(U1)  +AC  A  +  C(t) 


(III. 64) 


At  each  value  of  ij  ,  the  roots  of  the  above  equation  can  be  obtained  numeri¬ 
cally  using  a  root  search  technique.  Depending  on  the  value  of  \  ,  singular 
points  may  exist  inside  the  field  £  £  F  ,  With  the  correct 

value  of  A  satisfying  the  condition  f(0;1)  =  O  ,  the  roots  of  the  above 
equation  lie  outside  the  field  and  the  solution  is  free  from  singularities. 


CHAPTER  IV 


APPLICATIONS 

IV.  1.  ADIABATIC  POINT  EXPLOSIONS* 

i)  Introduction 

In  this  section,  a  detailed  solution  for  the  non-self-similar  blast  waves  in  an 
adiabatic  medium  is  obtained.  The  three  geometrical  symmetries  of  blast  waves: 
spherical,  cylindrical  and  planar  are  considered.  The  flow  field  is  considered 
as  a  sourceless  one,  while  the  medium  is  treated  as  a  perfect  gas  with  constant 
specific  heats.  The  solution  is  based  on  satisfying  the  condition  of  zero  par¬ 
ticle  velocity  at  the  center  of  symmetry  rather  than  the  constant  energy  crite¬ 
rion.  Using  this  condition,  a  considerable  saving  in  the  numerical  solution 
steps  is  achieved.  Finally,  the  accuracy  of  the  present  solution  is  discussed 
by  comparing  it  with  other  solutions. 

ii)  Problem  Formulation 

The  basic  equations  governing  any  gasdynamic  problem  are  given  by  Eqs.  (11.20), 
(11.22)  and  (11.27).  For  a  perfect  gas  with  constant  specific  heats,  a  source- 

less  flow  field  and  uniform  ambient  density,  these  equations  may  be  represented 
as: 

and 

iL)-M  =  o  J 

The  boundary  conditions  of  the  problem  are  given  by  the  Rankine-hugoniot  rela- 


*This  application  is  based  on  Attia  (1974), 
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tions,  Eqs.  (11.88)  -  (11.90),  which  may  be  rewritten  as  follows 


F„  =  f(b *)  =  (i  - />)  a  -  *) 


k  =  A-  (i,  tf)  s  - - — - 

p  +  a  -w 


and 

9„  =  3(1,3)  =  fl-P)  (1--JL  a) 


where 


/3  *  JLzl. 

avi 


The  application  of  the  quasi-similar  approximation  on  the  governing 
Eqs.  (IV. 1.1),  yields  the  following  ordinary  differential  equations 
by  Eqs.  (III. 23)  -  (III. 25)  when  all  <p  's  vanish 

(f-x)  ^  +  k(j£  +  i£)  +  \Wf\so 

<f-x>£+t£  +  UB-  °  r 

and 

(f.x)gL  +  X3(ML  +  jJL)  +  Mc  so  . 


(IV. 1.2) 


equations, 
as  given 


(IV. 1.3) 


where 
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A  =  CL^Ul. 

p  +(i~p)y 


0£. 


i/+l 


*(*-*) 


and 


c  «  - 


4  -  £  u/ 


Equations  (IV. 1.3)  may  be  put  in  the  following  autonomous  form,  as  given  by 
Eqs.  (III. 26)  -  (III. 28),  with  all  's  set  to  squal  zero 


and 


_  i|g£  +9(9  6  8 

k(f-xf-X3 

Ak  =  -zJl,  (i£  +  j.£  +  Ad) 

d*  <tx  * 


ii 

d* 


= .  /i  +  afs] 


(IV. 1.4) 


An  overall  check  for  the  denisty  profiles  is  the  mass  integral,  Eq.  (III. 65), 
which  can  be  rewritten  as 


=  / 
o 


k  xJ  <1  % 


s 


i 
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(IV. 1.5) 


Furthermore,  multiplying  the  first  of  Eqs.  (IV. 1.3)  by  and  arranging  terms, 
one  has 
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i+i 


dk 


i  % lmi  fk  +  akrt  xi  =0 


Integrating  the  above  equation  with  respect  to  X  from  X  =  0  to  X  =1,  one 
gets  after  some  algebraic  manipulations 


/  J(x‘fk)  .  /  d  (xulk)  +  (Ui  +  *X){  k  ix 


x  O 


which  gives 


f 

hn  (Pn„  i )  +  (kl+  AX)  /A  Jx  x  o 

O 


Using  Eqs.  (IV. 1.2),  jne  obtains 


i 

f  k  Xs  AX  x _ * _ 

°  Ui+AX 


(IV. 1.6) 


Consequently,  the  density  profile  obtained  by  the  quasi-similar  solution  does 
not  satisfy  the  global  mass  conservation.  However,  it  should  satisfy  the  quasi 
similar  mass  integral  given  by  Eq.  (IV. 1.6). 

The  difference,  £  ,  between  the  two  mass  integral  equations,  which  represents 

the  inherent  deviation  introduced  by  the  quasi-similar  approximation,  can  be 
obtained  from  Eqs.  (IV. I. 5)  and  (IV. 1.6)  as 


€  = 


4  +  i  +  A  %  4  •+  1 


or 


-AX 

V+i)(Ui  +  flX') 
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Since  both  ft  and  A  are  functions  of  ,  the  error  in  the  mass  integral  is 
also  a  function  of  y  .  It  is  obvious  that  the  error  is  zero  when  fj  =  0  and 
this  corresponds  to  =  0  which  is  the  self-similar  case. 

The  energy  integral,  Eq.  (11.95),  is  reduced  to 


] 


(IV. 1.7) 


while  the  decay  coefficient  A  takes  the  form  given  by  Eq.  (11.96),  which  is 


(Ui)  3 3  - 
j;  -  4 

di/ 


(IV. 1.8) 


THe  adiabatic  integral,  Eq.  (III. 3d),  which  is  used  as  a  check  on  the  precision 
of  the  numerical  solution,  may  be  rewritten  as  follows 


k!  [k  - i _ 

3  )  ffi  Ai-M] 

P  +1  *" 


(IV. 1.9) 


with 


^  ^  '£  -  if!) 

4+1  ♦ 


The  solution  is  considered  to  be  correct  when  it  satisfies  the  compatibility 
equation,  Eq.  (II. 11?)',  that  is 
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=0 


CIV. 1.10) 


iii)  Solution 

The  governing  equations  of  the  quasi-similar  solution  in  their  autonomous  form, 
Eqs.  (IV. 1.4),  which  are  coupled  ordinary  differential  ones,  are  subject  to  the 
boundary  conditions  given  by  Eqs.  (IV. 1.2).  .These  equations  can  be  easily  inte¬ 
grated  numerically,  once  the  correct  value  of  ^  at  a  specific  ij  is  determined. 
There  are  two  criteria  for  evaluating  the  shock  decay  coefficient  %  .  One  is 
based  on  the  conservation  of  total  energy  behind  the  shock  and  the  other  is 
based  on  satisfying  the  condition  of  symmetry,  namely,  f  =  0  at  the  center  of 
symmetry.  Although  only  the  latter  is  used  in  these  calculations,  it  may  be 
interesting  to  summarize  the  steps  of  calculations  which  one  may  follow  to  de¬ 
termine  the  correct  A  based  on  each  criterion. 

First  of  all,  one  should  chose  the  values  of  step  sizes  of  both  5C  and  3  .  The 
whole  range  of  3  (o-<y<i)  is  divided  into  small  intervals,  not  necessarily 
equal.  The  step  size  a  X  is  kept  constant  throughout  the  integration  inter¬ 
vals,  and  assigned  different  values  until  the  values  of  f  ,  (\  and  3  ,  ob¬ 
tained  from  the  numerical  calculations,  satisfied  the  adiabatic  integral,  Eq. 

(IV. 1.9),  with  the  desired  accuracy. 

If  the  solution  at  3*3;  is  known,  then  the  solution  at  = 

1l  +  C4«;  ),  based  on  the  constant  energy  criterion,  may  be  computed  as  fol¬ 
lows 


a)  At  the  value  of  3  =  3iV<  ,  the  boundary  conditions  f  (1,  3|Vi  ), 

k  (1,  *dj  £  )  and  9  (1,  )  can  be  computed  from  Eqs.  (IV. 1.2). 

b)  A  trial  value  of  the  decay  coefficient  \  is  then  chosen  and  Eqs.  (IV. 1.4) 
can  be  solved  numerically. 


c)  The  value  of  the  energy  integral  function  3^,  (  can  then  be  determined 
from  Eq.  (IV. 1.7)  using  any  numerical  integration  method,  for  example,  the 
Simpson's  rule. 


d)  ror  small  intervals  in  3 


the  derivative 


d 

d3 


in  Eq.  (IV. 1.6)  can  be 


put  in  a  finite  difference  form  as 
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B  A  ft  _  -  3j)i 


s  jgjllki  ~  ^0*' 
(A 


Thus  a  calculated  value  of  A  is  given  by 


(Jtl)  “  (V»Vl  /  itoVill 


Jj),']  /«  *>,■ 


e)  The  above  procedure  is  repeated  for  other  trial  values  of  X  and  a  curve 
between  the  calculated  values  and  the  assumed  ones  may  be  drawn.  The  correct 
value  of  X  will  be  obtained  as  the  intercept  of  the  45°  line  with  the  above 
curve  and,  with  this  correct  value,  the  profiles  of  the  various  fluid  proper¬ 
ties  are  then  determined. 


f)  Having  computed  everything  at  id  s  id,v<  *  onr>  can  proceed  further  until 
the  solution  is  obtained  for  the  whole  range  of  id  . 

I 

To  obtain  the  starting  conditions  for  this  computational  procedure,  one  simply 
substitues  id  "  0  at  the  system  of  Eqs.  (IV. 1.4).  They  become 


e i£  m 
JX 


k(f-x)1  ~  *9 


and 


.gLH.  -  -  A  f  JJL  +  j  j_ ) 
Jx  :f~x)  x  ' 


dl 


=  #  ] 


65 


Similarly,  the  substitution  of  if  =  0  into  Eq.  CIV. 1.8)  yields 


X=  3  +  1 


The  boundary  conditions. 


fn  - 


and 

%  -  ‘-t* 


Eqs.  (IV. 1.2),  become 


2 

nvi 


X+i 
2 r-i 


2 

V+i 


Fortunately,  this  system  of  equations  has  been  solved  exactly  using  the  concept 
of  self-similarity.  Thus,  the  value  of  (Tj)  (a  =  0)  can  be  determined  and 
is  used  as  the  starting  value  for  the  computational  procedure. 

On  the  other  hand,  the  solution  based  on  zero  particle  velocity  at  the  center 
of  symmetry,  which  is  applied  here,  required  the  following  calculation  steps: 

a)  At  any  value  of  ,  the  boundary  conditions  f(lsit)  ,  and 

4  a,  it)  can  be  computed  from  Eqs.  (IV. 1.2). 

b)  A  trial  value  of  the  decay  coefficient  A  is  assumed  and  Eqs.  (IV. 1.4) 
can  be  integrated  numerically .  Due  to  the  singularity  at  the  center  of  symme¬ 
try,  the  integration  is  stopped  at  X  r  X  which  is  arbitrarily  close  to  X  =  0. 

c)  The  velocity  at,  the  center  of  symmet:^  can  then  be  obtained  using  li¬ 
near  extrapolation  from  %  =  X  to  PC  =0. 

d)  1  he  above  procedure  is  repeated  for  several  values  of  assumed  decay  co¬ 
efficient  and  a  curve  is  plotted  between  the  assumed  value  of  A  and  the  velo- 
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city  at  the  center  of  symmetry.  The  correct  value  of  X  is  obtained  as  the 
intercept  of  the  above  curve  with  the  line  f  =  0. 

e)  With  the  correct  value  of  X  determined,  the  flow  field  parameters  are 
to  be  obtained  by  matching  the  numerical  solution,  as  X  ranges  over  from 
X  =  1  to  x  &  X  ,  with  appropriate  asymptotic  formulae  valid  from  X  r  9C 
to  X  =  0,  as  that  given  by  Eqs.  (III. 57]  -  (III. 59),  which  are 

f-£x 

k  =  kc  X* 

and 

9  =  9(0,10+  9„  X 

t 

The  value  of  fo  is  determined  as 

fc  C»)  . 

X 


consequently,  the  velocity  at  any  value  of  X  ^  X  is  given  by 


f  (X,  4)  =  .  EJ&li  .  x 

X 


une  may  now  proceed  to  uetermine  the  density  and  pressure  profiles  near  the 
uunter.  The  exponent  o£  can  be  determined  from  the  following  relation 
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hOj+li  -  (  #  f 

k  (x-  aXj3)  *  “ 

where  a#  is  the  step  size  used  in  the  numerical  iterat^v^  integration. 

Having  computed  o<  from  the  preceeding  equation,  the  value  of  the  function  h.c 
can  be  determined  from 

k.  (y)  s  JySlj.^1 

°  —  et 

CC 

Thus,  the  density  at  any  value  of  X  s?  X  can  be  obtained  as 

k  (%,3)  =  k0(H)  X* 


In  order  to  find  the  pressure  profile  near  the  center,  one  must  first  determine 
the  functions  S(Oj'd)  and  30C^)  •  They  can  be  obtained  as  the  solution 

of  the  two  simultaneous  equations 

3(Xj4)  =  S(o,<j)+  90C3)  x 

and 

9  (x  -  lx  ,  3)  =  3(0;  d)  +  y^cd)  (x-  ax)**2 

Solving  for  the  unknown  functions  $  (*J)  and  3(0)3)  ,  one  obtains 


3  (H)  =  9  (x-ax,3)  -  3(xjd) 

.r-  *+■  2  — .  oi>+  2- 

(X-aX)  -  X 
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and 


6(o,  = 


(x-ixf+i  SOejM-X**1  SCX'-AX,*') 


-  .ai+2  —  et+  2 

(X-AX)  -  X 


consequently  the  pressure  at  any  value  of  X  ^  X  is  given  by 


9  (XjV)  =  s (o,  <d)  +  6  ca;  # 


iv)  Results  and  Conclusions 

The  self-similar  solution  is  illustrated  in  Figs.  IV. 1.1  to  Fig.  IV. 1.4  j  Hcf4. 
It  is  seen  that  the  pressure  immediately  behind  the  shock  wave  is  a  maximum 
and  falls  off  quite  rapidly  near  the  shock  wave  to  a  nearly  constant  value  for 
%  <  0.5.  As  expected,  the  fall-off  in  pressure  is  greatest  for  the  spherical 
case,  with  its  greatest  freedom  for  expansion,  and  least  for  the  planar  case. 
The  same  tendency  is  even  more  accentuated  in  the  density  profiles,  where  it  is 
seen  that  nearly  all  of  the  mass  of  the  gas  engulfed  by  the  blast  is  concentra¬ 
ted  close  to  the  shock  front  itself.  The  effects  of  pressure  and  density  chan¬ 
ges  are  reflected  in  the  temperature  profiles,  since  T  «  P/j7  .  It  is  seen 
that,  subject  to  the  assumption  of  an  inviscid,  non-conducting  and  non-radia¬ 
ting  gas,  enormous  temperatures  are  developed  towards  the  center  of  the  blast, 
as  a  result  of  the  prevailing  vanishingly  small  densities  and  finite  pressure. 

I 

This  is  particularly  marked  in  the  spherical  case.  The  particle  velocities  de¬ 
crease  from  their  maximum  values  immediately  behind  the  shock  front  to  zero  at 
the  origin  of  the  blast.  It  is  observed  that  the  curves  differ  only  slightly 
in  the  range  0.5  <%<  1.  At  any  given  X  >  the  spherical  flow  velocity  is  the 
lowest  and  the  planar  velocity  the  largest. 

During  the  decay  of  the  shock  front  from  its  strong  limit,  self-similar  case, 
to  a  sound  wave  (  M  =  1),  the  quasi-similar  model  is  used  to  predict  the  non--- 


self-similar  case  based  on  satisfying  the  condition  of  zero  particle  velocity 
at  the  center.  This  criterion  is.  of  course,  not  new.  It  was  used  previous¬ 
ly  in  finding  higher  order  terms  in  the  perturbation  technique  (Korobeinikov 
and  Mel'nikova  (1962);  Bach  and  Lee  (1969)). 

Conditions  behind  the  shock,  as  given  by  the  solution  based  on  the  above  men¬ 
tioned  criterion,  are  shown  for  air  (  ■=  1.4)  in  Figs.  IV. 1.5  -  IV. 1.16: 

they  define  the  motion  and  the  changes  of  the  state  of  the  gas  in  the  three 
cases  for  J  =  0,  1,  2  and  for  different  values  of  the  shock  wave  Mach  number. 
It  is  observed  that  the  decrease  in  the  shock  Mach  number  causes  the  pressure 
to  be  uniform  in  the  whole  range  until  the  Mach  number  approaches  unity 
where  the  pressure  becomes  vanishingly  different  from  the  ambient  pressure. 

It  is  also  seen  that  at,  high  values  of  the  shock  Mach  number  (  fi  ^  3)  the 
density  falls  off  very  rapidly  behind  the  shock.  At  low  values  of  the  shock 
Mach  number,  the  mass  engulfed  by  the  front  becomes  equally  distributed  inside 
the  wave  and  the  density  approaches  the  ambient  density  everywhere  except 
at  X  -  0.  Particle  velocities  decrease  as  the  time  passes  and  tend  to  zero 
when  the  shock  wave  attenuates  to  a  sound  wave. 

The  usefulness  of  the  derived  asymptotic  formulae,  Eqs.  (III. 57)  -  (III. 59), 
is  apparent  in  the  neighbourhood  of  X  ~  0,  where  the  original  equation  shows 
a  singularity.  In  this  region,  the  velocity  distribution  is  linear.  Close  to 
zero,  the  pressure  is  finite  and  asymptotically  constant  to  the  X  coordinate. 
The  density  approaches  zero  very  rapidly  indicating  that  the  gas  is  displaced 
from  the  center  of  the  explosion.  The  temperature  distribution  T p//  shows 
a  steep  increase  to  infinite  values. 

In  the  course  of  high  temperature  processes  within  the  blast  the  specific  heat 
ratio,  ,  differs  from  its  initial  value  of  1.4.  Therefore,  the  solution  is 
also  obtained,  when  assuming  air  to  behave  as  a  perfect  gas  with  a  constant 
average  value  of  specific  heat  ratios  of  1.2  and  1.3.  Figures  IV.1.17-IV.1.20 
:how  a  comparison  between  the  distribution  of  fluid  properties  behind  the  wave 
front  for  the  special  case  of  J  =  1,  H  =5.  The  trend  is  the  same  for  ether 
values  of  fl  as  well  as  for  other  types  of  blast  waves  (i.e.  i  =  0  and  2). 

The  accuracy  of  satisfying  the  adiabatic  integral,  Eq.  (IV. 1.9)  is  the  crite¬ 
rion  for  the  choice  of  the  step  size  a  X  in  the  Runge-Kutta  iterative  inte- 
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gration.  Figure  IV.  1.21  shows  that  for  =  0.01,  the  adiabatic  integral  is 

satisfied  within  1  x  10  It  is  also  seen  that  at  first  the  error  decreases 
quite  rapidly  with  the  step  size  and  then  becomes  almost  constant  for 
0.005. 

&3  C 

The  integral  curves  in  Z  -  F  phase  plane  where  Z  s  and  f  e.  X  ,  are 

shown  in  Figs.  (IV.1.22)(  -  (IV.  1.24)  for  different  values  cf  the  shock  Mach 
number.  All  the  curves  start  from  the  Rankine-Hugoniot  curve  and  all  of  them 
tend  to  infinity  corresponding  to  ^  =  0.  It  is  seen  that  the  locii  of  the 
singularity  2  =  (1  -  r  )  ,  Eq.  (III. 61),  and  the  singularity  F  =  1,  Eq. 

(III. 00),  do  not  intersect  with  any  of  the  integral  curves,  i.e.  the  solutions 
are  singularity  free.  Of  course,  all  the  integral  curves  are  drawn  with  the 
correct  value  of  X  satisfying  the  condition  ft*,  tf) -  0. 

The  effect  of  inserting  values  of  ^  in  the  differential  equation  that  are 
less  or  greater  than  the  correct  one  is  indicated  in  Fig.  (IV. 1.25)  where  it 
is  seen  that  for  'X  less  than  the  correct  one,  the  singularity  f  -  1  is  en- 

a  2 

countered  while  for  A  greater  than  the  correct  value,  the  singularity  2  s=  (i-F) 
is  reached. 

It  may  be  observed  from  Fig.  IV. 1.26  that  a  considerable  error  is  obtained  in 
the  mass  integral  which  means  that  the  density  profiles  are  not  accurate. 
However,  the  comparison  made  in  Fig.  IV. 1.27  between  the  present  solution  and 
the  solution  of  Bach  and  Lee  (1969)  indicates  that  the  present  solution  has  a 
wider  range  of  applicability  and,  therefore,  it  gives  a  reasonable  description 
of  the  shock  trajectory. 

The  Rankine-Hugoniot  shock  conditions,  which  represent  the  maximum  damage  oc¬ 
curring  due  to  the  passage  of  tne  blast  wave,  are  shown  in  Figs.  IV. 1.20  - 
IV. 1.31.  It  is  observed  that  when  the  blast  expands,  the  pressure,  tempera¬ 
ture  and  density  behind  the  shock  wave  approach  the  atmospheric  conditions 
while  the  velocity  behind  the  shock  wave  tends  to  zero.  The  relatiun  of  the 
shock  wave  attenuation  is  shown  in  Fig.  IV. 1.32,  where  it  is  seen  that  the 
shock  Mach  number  approaches  unity  when  the  shock  radius  tends  to  infinity. 

The  relation  between  the  shock  radius  and  the  time  is  shown  in  Fig.  IV. 1.33. 

The  decay  coefficient  obtained  by  the  present  criterion  is  shown  in  Fig. 

IV. 1.34  as  a  function  of  ^  for  the  three  types  of  blast  waves:  spherical,  cy- 
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lindrical  and  planar. 

The  present  calculation  is  compared  with  that  of  Dshima  (1964)  in  the  A  -  ^ 
curve,  as  shown  in  Fig.  IV. 1.35.  It  is  seen  that  the  present  results  give  a 
slightly  greater  value  for  the  decay  coefficient.  However,  since  the  flow 
field  is  very  sensitive  to  A  >  this  slight  change  may  change  the  flow  field 
considerably. 

It  should  be  noted  that  either  the  zero  particle  velocity  criterion,  or  the 
energy  integral  criterion  is  sufficient  in  iterating  for  i.he  correct  solution. 
However,  it  was  found  that  the  solution  based  on  zero  particle  velocity  at  the 
center  is  by  far  more  sensitive  than  the  energy  integral  to  very  small  varia¬ 
tions  in  A  as  it  approaches  its  correct  value.  Hence,  highe.  accuracies  in 
the  solutions  can  be  obtained  using  the  zero  particle  velocity  criterion.  Fi¬ 
gure  IV. 1.36  shows  the  variation  of  the  particle  velocity  at  the  center  fa,*) 
with  the  assumed  value  of  A  •  It  is  observed  that  any  small  change  in  the 
value  of  the  decay  coefficient  produces  a  considerable  change  in  the  velocity 
at  the  center.  In  the  present  calculations  the  velocity  at  the  center  is  con¬ 
sidered  zero  when  ]FeM  «  io'4. 

The  iteration  by  the  energy  integral  criterion  is  not  rapidly  convergent.  Fur¬ 
thermore,  when  A  approaches  its  correct  values,  the  properties  of  the  solu¬ 
tion  by  the  zero  particle  velocity  criterion,  such  as  the  position  of  the  sin¬ 
gular  points  and  whether  the  velocity  at  the  center  is  positive  or  negative, 
may  indicate  the  direction  for  the  correct  values  of  the  decay  coefficient.  In 
addition,  reducing  the  step  of  the  numerical  integration,  results  in  a  consi¬ 
derable  saving  in  machine  time. 

The  difference  between  the  results  obtained  from  the  two  criteria  is  due  to 
the  quasi-similar  approximation  which  reduced  the  partial  differential  equa¬ 
tions  to  oroinary  differential  ones.  It  should  be  noted  that  if  the  governing 
equations  were  solved  exactly,  the  two  criteria  would  result  in  the  ^ din> ■■  va¬ 
lues  of  the  decay  coefficient  and  thu  flow  field  variables. 


Fig.  IV. 1.1  Self-Similar  pressure  profiles  for  intense  planar,  cylindrical 
and  spherical  blast  waves  with  IS  =  1.4. 

Fig.  IV. 1.2  Self-similar  density  profiles' for  intense  planar,  cylindrical 
and  spherical  blast  waves  with  =  1.4. 

Fig.  IV. 1.3  Self-similar  temperature  profiles  for  intense  planar,  cylindri¬ 
cal  and  spherical  blast  waves  with  if  =  1.4. 

Fig.  IV. 1.4  Self-similar  velocity  profiles  for  intense  planar,  cylindrical 
and  spherical  blast  waves  with  -  1.4. 

Fig.  IV. 1.5  Pressure  distribution  behind  the  wave  front  for  spherical  blast 
waves  at  different  values  of  the  shock  Mach  number  hi  while 

*  =  1.4. 

Fig.  IV. 1.6  Density  distribution  behind  the  wave  front  for  spherical  blast 
waves  at  different  values  of  the  shock  Mach  number  m  while 
=1.4.’ 

Fig,  IV. 1.7  Temperature  distribution  behind  the  wave  front  for  spherical 
blast  waves  at  different  values  of  the  shock  Mach  number  M 
while  =  1.4. 

Fig.  IV. 1.8  Particle  velocity  profile  behind  the  wave  front  for  spherical 
blast  waves  at  different  values  of  the  shock  Mach  number  ^ 
while  )S  =  1-4. 

Fig.  IV. 1.9  Pressure  distribution  behind  the  wave  front  for  cylindrical 
blast  waves  at  different  values  of  the  shock  Mach  numbe  I'  M 
while  Y  =  1.4. 
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Fig.  IV. 1.10  Densicy  distribution  behind  the  wave  front  for  cylindrical 
blast  waves  at  different  values  of  the  shock  Mach  number  Ml 
while  tf  =  1.4. 

Fig.  IV. 1.11  Temperature  distribution  behind  the  wave  front  for  cylindri¬ 
cal  blast  waves  at  different  values  of  the  shock  Mach  number 
M  while  8"  =  1.4. 

Fig.  IV. 1.12  Particle  velocity  profile  behind  the  wave  front  for  cylindri¬ 
cal  blast  waves  at  different  values  of  the  shock  Mach  number 
Ml  while  If  =  1.4. 

Fig.  IV. 1.13  Variation  of  the  pressure  profile  with  the  shock  Mach  number 
Ml  for  planar  blast  waves  while  "S'  -  1.4. 

i 

Fig.  IV. 1.14  Variation  of  the  density  profile  with  the  shock  Mach  number 
t-|  for  planar  blast  waves  while  Jf  =  1.4. 

Fig.  IV. 1.15  Variation  of  the  temperature  profile  with  the  shock  Mach  num¬ 
ber  for  planar  blast  waves  while  S'  =  1.4. 

Fig.  IV. 1.16  Variation  of  the  particle  velocity  profile  with  the  shock  Mach 
number  ft  for  planar  blast  waves  while  S'  =  1.4. 

Fig.  IV. 1.17  Comparison  between  the  pressure  profiles  for  different  values 
of  the  specific  heat  ratio  while  J  =  1  and  Ml  =  5. 

Fig.  IV. 1.10  Comparison  between  the  density  profiles  for  different  values 
of  the  specific  heat  ratio  while  J  =  1  and  Ml  =  5. 

Fig.  IV. 1.19  Comparison  between  the  temperature  profiles  for  different  va¬ 
lues  of  the  specific  heat  ratio  while  J  =  1  and  M  =  5. 
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Fig.  IV. 1.20  Comparison  between  the  particle  velocity  profiles  for  diffe¬ 
rent  values  of  the  specific  heat  ratio  while  J  =  1  and 

!“1  =  5. 


Fig.  IV. 1.21  Error  in  the  machine  computations  as  a  function  of  the  step 
size  in  the  Runge-Kutta  iterative  integration  in  the 

case  of  J  =  1,  H  =  2  and  -  0.1. 


Fig.  IV. 1.22  Effect  of  variation  of  the  shock  Mach  number  M  on  the  loca¬ 
tion  of  the  integral  curves  with  respect  to  the  locii  of  the 

2 

singularities  f  =  1  and  Z  =  (1  -  F  )  for  d  =  2  and 
=  1.4. 


Fig.  IV. 1.23  Effect  of1  variation  of  the  shock  Mach  number  M  on  the  loca¬ 
tion  of  the  integral  curves  with  respect  to  the  locii  of  the 

r-  2 

singularities  f  =  1  and  Z  =  (1  -  r  )  for  d  =  1  and 
=1.4. 


Fig.  IV. 1.24  Effect  of  variation  of  the  shock  Mach  number  M  on  the  loca¬ 
tion  of  the  integral  curves  with  respect  to  the  locii  of  the 
singularities  F  -  1  and  Z  =  (1  -  F  )  for  d  =  0  and 
=1.4. 

2 

Fig.  IV. 1.25  Location  of  the  singular  points  F  =  1  and  Z  =  (1  -  F  )  at 
values  of  rather  than  the  exact  one  in  the  case  of  d  -2 
and  M  =  2. 


Fig.  IV. 1.26  The  mass  integral  (  J  +  1)  as  a  function  of  in  the 
case  of  If  =  1.4. 

Fig.  IV. 1.27  Comparison  between  the  present  solution  with  the  first  order 
perturdation  solution  (Sakurai,  1954)  and  the  second  pertur¬ 
bation  solution  (Bach  and  Lee,  1969). 
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Fig.  IV. 1.28  Pressure  ratio  across  the  shock  front  as  a  function  of  the 
non-dimensional  shock  radius  J  for  the  three  cases  of  sphe¬ 
rical,  cylindrical  and  planar  blasts  while  X  -  1.4. 

Fig.  IV. 1.29  Density  ratio  across  the  shock  front  as  a  function  of  the 

non-dimensional  shock  radius  J  for  the  three  cases  of  sphe¬ 
rical,  cylindrical  and  planar  blasts  while  5  =  1.4. 

Fig.  IV. 1.30  Temperature  ratio  across  the  shock  front  as  a  function  of 
the  non-dimensional  shock  radius  J  for  the  three  cases  of 
spherical,  cylindrical  and  planar  blasts  while  =  1.4. 

Fig.  IV. 1.31  Variation  of  the  particle  velocity  at  the  shock  front  as  a 
function  of  the  non-dimensional  shock  wave  radius  J  for 
the  three  cases  of  spherical,  cylindrical  and  planar  blasts 
while  =  1.4. 

Fig.  IV. 1.32  Mach  number  of  spherical,  cylindrical  and  planar  blast  wave 
front  as  a  function  of  the  non-dimensional  radius  while 
=  1.4. 

Fig.  IV. 1.33  Variation  of  the  shock  wave  non-dimensional  radius  J  with 
the  non-dimensional  time  for  spherical,  cylindri¬ 

cal  and  planar  blasts  with  if  =  1.4. 

Fig.  IV. 1.34  Shock  decay  coefficient  ^  as  a  function  of  Jj  for  the  three 
cases  of  spherical,  cylindrical  and  planar  blasts  with 

#  =1.4. 

fig.  IV. 1.35  Comparison  between  the  solution  based  on  the  present  crite¬ 
rion  with  the  solution  based  on  the  energy  integral  criterion 
for  j  -  1  and  £  =  1.4. 

Fig.  IV. 1.36  Variation  of  the  velocity  at  the  center  with  the  assumed  ^alue 
of  the  shock  decay  coefficient  in  the  case  of  J  -  1  and 

H  =  2. 
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IV.  2.  QUASI-SIMILAR  solutions  of  adiabatic  point  explosions— 

AN  ANALYTICAL  SOLUTION*  7 

i)  Introduction 

In  this  section,  the  quasi-similarity  is  introduced  to  study  the  flow  field 
structure  of  adiabatic  blast  waves.  The  formulation  of  the  problem  is  made 
applicable  for  the  three  geometrical  symmetries:  spherical,  cylindrical  and 
planar  waves.  The  flow  field  is  assumed  to  be  inviscid  and  the  flowing  me¬ 
dium  is  treated  as  a  perfect  g'.s  with  a  mean  specific  heat  ratio  .  An  ap¬ 
proximate  analytical  solution  is  obtained  and  the  results  are  compared  with 
other  existing  solutions.  Finally,  depending  on  these  results,  the  validity 
af  th8  quasi-similar  theory  is  discussed. 

i i )  Problem  Formulation: 

The  conservation  equations  in  their  general  form  are  given  by  Fqs.  (11.20), 
(11,22)  anu  (11.24).  /or  inviscid  and  sourceless  flow  field,  all  <p>t  in 
these  equations  vanish.  Moreover,  when  the  flowing  medium  js  treated  as  a 
perfect  gas  with  a  mean  specific  heat  ratio,  ,  and  when  P  is  constant, 
these  equations  are  simply  reduced  to 


J  2k  *  ~*k.  +  knJL  t  il)  =r o 

x  v  -v/e  x 


*This  application  is  paced  cn  Osfima  (1964). 
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The  governing  equations,  Eqs.  (IV. 2.1),  are  subject  to  the  boundary  conditions 
given  by  the  Rankine-Hugoniot  relations,  Eqs.  (11.88)  -  (11.90),  which  may  be 
rewritten  as  follows 


fa  =  PQ'T>  =  d- w 

0  +  1 


9„  =  9  d,  J  >  =  .KPr^zilM. 


(IV. 2. 2) 


H„  =  kd,  J)  = - **  i  — ■ 

y-l+d 


For  adiabatic  point  explosions  and  sourceless  flow  fields,  the  energy  integral, 
Eq.  (11.95),  reduces  to 


_  -  (J+  1) 

<T  -  JL  T  7  -e _ t _ -  1 

3  y  *■  J  (/+.<>  (»•- a;  j 


(IV. 2. 3) 


Utilizing  the  concept  of  quasi-similarity,  the  governing  equations  reduce  to  a 
system  of  ordinary  differential  ones,  given  by  Eqs.  (III. 23)  -  (III. 25)  in 
which  all  the  ’ s  equal  to  zero  for  a  sourceless  f 3 field.  Therefore,  one 
may  rewrite  the  governing  equations  of  the  quasi-sir  liar  solutions  as  follows: 


(f.x)  4£  +.L  dl  +  ief  =o  ,  uv.2.4) 

cJ  /t  Or* 


(f-x)4 -  +  -° 
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where 


2d 

i  +  2  y 
^  +  1 

2  (i-d) 


a-x  -  <ir-d)  y 


The  adiabatic  integral,  Eq.  (III. 33),  may  be  rewritten  directly  as 

k  [k(x-f)%6]  J  +  t  -***/(*-i  +  t*) 

.  ycy+pQ+OA' y-  i+2to] 

Ay-  in/  + 


(IV. 2. 5) 


In  order  to  use  this  integral  as  a  check  through  the  numerical  integration, 
one  may  introduce  the  following  new  variables 


I  =  *  -  f 


K  «r  1  - 

ys 


(IV. 2. 6) 


then,  the  basic  equations,  Eqs.  (IV. 2. 4),  in  terms  of  I  ana  K  ,  will  take 
the  form 
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dJi  *  l=A  f  *ji£.AI~  x(i-M-)  -  —  -^-1 

d%  i  l  3  dx  K  7*'  k  Jjc  J 

.  ij* 


where 


i.iiilL.  A 

*<*-*) 

q  -  r  j_±j* _ £ _ 1  ^  +  i 


C  c  LL±JQl.  ^ 


and 


D/  -  1)  (iVi)  -h  -  ■ 1  ^  — — — 

#- 1+  *3 


1*1 

X  €  -  (X- 1)  ij 


and  the  boundary  conditions  are  transformed  to 


X„  *  1(1,1)  »  i-  t  =  -£0±li 

n  *  8>i 


Kn  =  K  (A/  *) 


«  j£±i>IL^ 

2  Jf_  (y.i)  y 


y  (IV. 2. 8) 


Thus,  the  problem  is  transformed  to  a  boundary  value  one  in  which,  A  is  an  un¬ 
known  factor  for  each  value  of  'd  .  An  iteration  procedure  is  performed  until 
the  correct  X  la  obtained,  which  must  satisfy  the  adiabatic  integral,  Eq. 
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(IV. 2. 5),  to  a  sufficient  accuracy.  The  values  of  X  determined  are  presented 
in  Table  IV. 2.1 

iil)  Approximate  Relations 

It  is  found  from  the  calculated  results  that  the  solution  may  be  represented 
approximately  in  the  form 

Is  X  +  p  Xn  (IV.  2. 9) 


where  oi  ,  (3  and  n  are  determined  by  the  boundary  conditions  and  the  asymp¬ 
totic  behaviour  of  I  and  K  to  X  -  0.  From  Eq.  (IV. 2. 7),  one  has 


X  •*  o 


L  im 

■x—o 


Lim  f A*+ 

X*~o 


thus 


iii  = 

<iX  X  —  O 


while  from  Eq.  (IV. 2. 9) 


|  =  <* 

X*o 


1  hus 


o( 


£+  e' 

f+i 


=  i  -  (J±—  )  ( — — o. - — 

v  J+l  '  28*  +  a  -  *3 


i  IV.  f.  I'll 
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Also,  from  Eq.  (IV. 2. 9),  one  has 


a  +  P 


Thus 


B  -  1  cx  -  y- 1  +  1 1  + _ — _  -  1  'IV. 2. 11) 

n  x+1  (j+i) 


while 


(IV. 2. 12) 


where 


(il )  s 

'ft 


_  il 

Jx 


✓ 


By  using  the  values  of  X  determined  previously,  the  values  of  ci  ,  and 
rt  can  be  obtained  from  the  above  relations.  These  values  as  well  as  tXj 
and  J  are  tabjlated  in  Table  IV. 2.2. 

"or  integrating  the  second  of  Eqs.  (IV. 2. 7),  it  follows  that 


i~  K 


JL  d  *  -  ShlL 

i  i 


d  x  -  J-&J2  Jx 
oc 


integration  of  the  above  equation  ;i,J,'. 


,  ,  ,  .  JT+  1  t 

,  /n  (i~  K)  =  0  / _ - -  /n  I  _  # 

V  «?(  +  9  o;R 


j  Or-  V 


PC 


„  rt- 2  „  ,  r  /*■!  i(Jf~  O  -j 

f  f  J_  „  r*L _ _ — .  1  Jit  -  f  n  [J  ^  j 

^  *■  *  (rt  -  i  J  (a:  +  p  n" 1 }  J 
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if*  1  /(*-*) 


-  fn  ( Cjnst-i ) 


Co  ft  ft J 


** 1  ♦  5^b&  t  %J(*'°~  *(n~x> 


The  constant  of  integration  is  determined  by  using  the  boundary  conditions  at 
^  =  1  to  yield  finally 


*n 


*VJ  ♦ 


o  •  no 

•Ufi-i)  -  *  /n_ 


fn-i) 


(IV. 2. 13) 


From  Eqs.  (IV. 2. 4),  one  has 


-li.f.  — ]J« 


which  yields  upon  integration 


J+i  - 


i  - - — X?  *  t  *  +  l 

k  =  Const.  ^  (X) 


n  (J+i  -  y-i  +  *!W  _  j 


Applyir<>  the  boundary  conditions  yields 


=  K  (-f) 


tM 

if- i  tid 
tx(r)-i) 


"(^-■gr T7i?)  ; 


(IV. 2.  l-i) 
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It  follows  also  from  Eqs.  (IV.2.4)  that 


Equations  (IV. 2. 9),  (IV. 2. 14)  and  (IV. 2. 15)  are  the  approximate  relations  for 
deternining  the  non-dimensional  gasdynamic  parameters  f  ,  k  and  6  ,  respec¬ 
tively. 

In  order  to  determine  a  relation  between  A  and  ,  one  may  differentiate 
the  relation  of  Eq.  (IV. 2. 9)  twice  and  eliminate  p  to  obtain  tt  -  1: 


(IV. 2. 16) 
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the  aid  of  Eqs .  CIV.2.BJ  and  CIV. 2. 9),  one  has 


* 


y-  i  2.  ‘d 


(  ~  <LzJ2 


(-@1- 

+  (Krt“d)  f1-  } 

I 

(ii.)  s  iA - [  <7-  ■■  <  <*-<>*«] 


.  ,  /  dLOL)  (-*!•£-)  »  X„  and  int0  Ec1, 

Substituting  the  values  of  \  ^X1  *  •  '  dX  n 

(IV.  .16).  a  third  order  algebraic  equation  In  X  It  obtained.  In  the  case 
of  p  anar  waves,  this  equation  Is  reduced  to  a  second  order  algebraic  equa- 
tion  and  its  proper  root  is  given  by 


X  -  - £r. 4  ** 

z.  <p 


(IV. 2. 17) 


where 


<p  =  D  F  j  V/=CF+OF+Q  j  0  ~  C  F  ~  A 

A*=  A  L4V-2  +  <V+5)y-(X-i)y2]  (l-'jfcr-i  +  Z'd) 

q" -  [4^-^  +  y  -  Of-0  tf2]  (*- *  +  *J0 

£•"=  -2  (1-  u)  px  -  (*-*)*] 

0^  s  _  4  (*V  i) 

r"=  -4  [>*-*  +  y  -  (*- 1)  y* ]  +  x  c*+  o  ^ 

r“-  *  v-m.?.**- -77f^7j  +  ^]-e-«^> 


iv)  Results  and  Conclusions: 

The  distributions  of  the  gasdynamic  parameters  £  ,  3  and  K  are  presented 
in  Figs.  IV. 2.1  to  IV. 2. 6  to  show  the  applicability  of  the  approximate  formu¬ 
lae,  given  by  Eqs.  (IV. 2. 9),  (IV. 2. 143  and  (IV. 2. 15),  which  are  compared 
with  numerical  and  Sakurai’s  second  approximation  perturbation  solutions. 
These  are  for  the  case  of  cylindrical  waves,  as  an  example,  with  two  diffe 
rent  values  of  Mach  number,  M  =  2  and  M  =  3,  while  X  =  1.4.  This  compa 
rison  proves  the  usefulness  of  these  formulae.  Especially,  in  the  neighbor 
hood  Cif  X  -  0,  where  the  original  equations  show  a  singularity  at  the  center 
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of  symmetry,  X  =  0,  these  formulae  yield  good  insight  into  the  problem  under 
study. 

The  pressure  and  temperature  distributions  are  also  shown  in  Figs.  IV. 2. 7  to 
IV. 2. 14,  for  the  three  geometrical  symmetries.  These  figures  show  that  the 
pressure  is  approximately  uniform  near  the  center  of  symmetry  while  the  tem¬ 
perature  has  a  steep  increase  near  the  center  of  symmetry  to  infinity  as  a 
result  of  neglecting  the  heat  transfer  effects. 

Figure  IV.2.15  shows  the  variation  of  X  with  for  different  values  of  ^  , 

in  the  case  of  planar  waves,  by  using  the  approximate  analytical  formula  gi¬ 
ven  by  Eq.  (IV. 2. 17). 

It  is  observed  from  the  experimental  data  (Oshima,  1960)  that  the  quasi-simi¬ 
larity  assumptions  fail  for  the  strong  blast  waves.  It  is  also  clear  that  the 
central  region  in  the  quasi-similar  solution  has  an  infinitely  high  tempera¬ 
ture  core  which  is  never  realized  in  the  actual  flow.  This  last  phenomena  is 
due  to  the  fact  that  transport  processes  have  been  neglected  in  the  governing 
Eqs.  (IV. 2.1)  and  not  due  to  quasi-similarity  assumptions.  Later  in  this 
chapter  we  shall  deal  with  this  issue. 

Therefore,  the  region  in  which  the  quasi-similar  solution  is  applicable  is 
limited  to  the  flow  fields  of  blast  waves  with  moderate  strength  excluding  the 
central  region. 

Figure  IV. 2. 16  shows  a  comparison  of  the  values  of  X  obtained  by  the  quasi¬ 
similar  technique  with  that  calculated  previously  by  SaKurai  (1961),  Sedov 
(195/)  and  Brode  (1955,  1959),  all  for  the  spherical  waves  with  i  =  '2. 

Sedov  (1957)  expressed  the  velocity,  pressure  and  density  as  functions  of  % 
and  b  .  These  variables  were  expanded  in  a  power  series  of  ,  and  J  was 
assumed  to  taKe  t hie  form. 


then  neglecting  the  higher  order  teims,  the  basic  equations  were  reduced  to  u 
system  of  linear  differential  equations  containing  the  unknown  constant  f\  . 
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This  value  of  /?  was  determined  so  that  the  boundary  conditions  at  X  -  1  and 
X  =  0  were  satisfied.  It  was  close  to  1.92  in  the  case  of  1=2.  Finally, 
the  results  obtained  were  expressed  by 

1  -  1 

1+1-922 


-  l.ss?  2  e 


1.92  2 


Brode  (1955,  1959]  calculated  numerically  the  adiabatic  point  explosion  problem 
using  artificial  viscosity  technique,  and  presented  an  approximate  formula  well 
in  agreement  with  the  calculated  results  using  the  quasi-similar  technique, 
which  is  given  by 


(1  +  2) 


f  3) — 2_ 


Finally,  Sakurai  (1954]  obtained  the  following  results 


A  =  3  (1  -  1-219  2) 


J3  =(!.<?#)  - - r— 

(1  -  1.918  2) 


which  is  also  compatible  with  previous  ones. 


< 
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Figure  Captions 

Fig.  IV. 2.1  The  distribution  of  the  non-dimensional  particle  velocity 

f  based  on  the  present  approximate  formulae  as  compared 
with  those  of  numerical  ahd  Sakurai's  second  order  pertur¬ 
bation  solutions  for  the  case  of  J  =  1  at  H  =3. 

Fig.  IV. 2. 2  The  distribution  of  the  non-dimensional  pressure  3  based 
on  the  present  approximate  formulae  as  compared  with 
those  of  numerical  and  Sakurai's  second  order  perturba¬ 
tion  solutions  for  the  case  of  J  =  1  at  *i  =  3. 

i 

Fig.  IV. 2. 3  The  distribution  of  the  non-dimensional  density  h  based 
on  the  present  approximate  formulae  as  compared  with 
those  of  numerical  and  Sakurai's  second  order  perturba¬ 
tion  solutions  for  the  case  of  J  =  1  at  M  =  3. 

Fig.  IV. 2. 4  The  distribution  of  the  non-dimensional  particle  velocity 
f  based  on  the  present  approximate  formulae  as  compared 
with  those  of  numerical  and  Sakurai’s  second  order  pertur¬ 
bation  solutions  for  the  case  of  J  *  1  at  fl  =  2. 


Fig.  IV. 2. 5  The  distribution  of  the  non-dimensional  pressure  3  based 
on  the  present  approximate  formulae  as  compared  with 
those  of  numerical  and  Sakurai’s  second  order  perturba¬ 
tion  solutions  for  the  case  of  i  =  1  at  H  =  2. 


Fig.  IV.2.6  The  distribution  of  the  non-dimensional  density  k  based 
ori  the  present  approximate  formulae  as  compared  ith 
those  of  numerical  and  Sakurai's  second  order  per.urba- 
tion  solutions  for  the  case  of  J  =  1  at  M  =  2. 


Fig.  IV. 2. 7  Pressure  distribution  at  different  values  of  the  shock 
Mach  number  fi  for  planar  waves  with  &  =  1.4. 


Fig.  IV. 2.  8  Pressure  distribution  at  ditr°rent  values  of  the  shock 
Mach  number  M  for  planar  waves  with  X  =  1.667. 

Fig.  IV. 2.  9  Pressure  distribution  at  different  values  of  the  shock 
Mach  number  M  for  cylindrical  waves  with  X  =1.4. 

Fig.  IV. 2. 10  Pressure  distribution  at  different  values  of  the  shock 
Mach  number  Ml  for  spherical  waves  with  X  =  1.4. 

Fig.  IV. 2. 11  Temperature  distribution  at  different  values  of  the 
shock  Mach  number  H  for  planar  waves  with  X  =  1.4. 

Fig.  IV. 2. 12  Temperature  distribution  at  different  values  of  the 

shock  Mach  number  Ml  for  planar  waves  with  "X  =  1.667. 

Fig.  IV. 2. 13  Temperature  distribution  at  different  values  of  the 

shock  Mach  number  Ml  for  cylindrical  waves  with  X  =1.4. 

i 

Fig.  IV. 2. 14  Temperature  distribution  at  different  values  of  the 

shock  Mach  number  Ml  for  spherical  waves  with  X  =  1.4. 

Fig.  IV. 2. 15  Variation  of  the  decay  coefficient  with  *d  in  the 
case  of  planar  waves  while  X  ■  1.2,  1.4  and  1,667. 

Fig.  IV. 2. 16  Comparison  of  the  values  of  \  based  on  the  present  ap¬ 
proximate  solution  with  those  of  Sakurai's  (1954), 

Sedov's  (1957)  and  Brode’s  (1955,  1959)  for  the  spheri¬ 
cal  waves. 
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IV.3  BLAST  WAVES  IN  REAL  GASES* 

i)  Introduction 

In  the  prevailing  number  of  cases,  with  very  few  exceptions,  despite  the  high 
temperatures  involved,  the  gaseous  medium  affected  by  blast  waves  has  been 
treated  in  the  literature  as  a  perfect  gas  with  constant  specific  heats.  How¬ 
ever,  due  to  the  extremely  high  temperatures  prevailing  close  to  the  centers 
of  blast  waves,  one  should  expect  the  excitation  of  all  possible  degrees  of 
freedom  of  the  internal  energy  of  the  gas.  One,  therefore,  should  take  into 
account  such  phenomena  as  vibrational  excitation,  dissociation,  electronic 
exitation  and  ionization,  as  well  as  the  influence  of  compressibility  at  high 
pressures. 

In  this  section,  the  real  gas  effects  in  adiabatic  point  explosions  are  taken 
into  account.  For  this  purpose,  an  equilibrium  thermodynamic  analysis  is  car¬ 
ried  out  in  different  ranges  of  pressures  and  temperatures.  At  nigh  pressures, 
the  compressibility  effect  is  considered  using  empirical  formulae,  such  as  the 
Beattie-Bridgsman  equation,  as  well  as  the  concept  of  residual  properties.  At 
high  temperatures,  on  the  other  hand,  quantum  statistical  thermodynamic  con¬ 
cepts  are  used  to  determine  the  degree  of  excitation  leading  to  dissociation 
and  ionization  of  the  gas.  As  formulated  previously  in  Chapter  II,  the  real 
gas  characteristics  are  expressed  in  terms  of  two  thermodynamic  properties. 
These  are  a  non-dimensional  speed  of  sound  factor,  /”"*  ,  and  the  dimensionless 
internal  energy,  s'  .  The  geometry  of  the  generated  field,  whether  it  is  a 
plane,  cylindrical  or  spherical,  is  considered.  The  flow  field  is  assumed  to 
be  sourceless,  inviscid,  non-conducting  and  non-radiating. 

ii)  Problem  Formulation 

The  general  form  of  the  governing  equations,  including  real  gas  effects,  are 
given  oy  Eqs.  (11.20),  (11.22)  and  (11.27).  In  the  absence  of  all  source 
terms  and  when  is  constant,  these  equations  become 


*  This  application  is  cased  on  Ghonier  (1975). 
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and 


.  J_  ( 


~>£if  +  $  -a  It  9 

kxf  ~*fcX’ 


-  o 


+ 


-*Lx, 


(IV. 3.1) 


where  f*  is  defined  by  Eq.  (II. 6)  as 

r  m  (ikl  \  f  a} 

1  (^T?  \  a  T  (IV*3 

The  boundary  conditions  of  the  problem  ax'e  given  by  Eqs.  (11.54)  -  (11.56) 
which  may  be  rewritten  as  fallows 


a  - L 


9  -  P  +  JL 

n  —  +  -=■ 


n  r 

lCL 


and 


**  S<sa  +■  -FT*  £  *  f  £ 
'a 


(IV. 3. 3) 


where  ok  is  defined  in  Eq.  (11.15)  as 


146 


C'  S  — 
Wn 


(IV. 3. 4) 


The  mass  and  energy  integrals,  Eqs.  (11.65)  and  (11.73),  in  this  case,  are  gi¬ 
ven  by 


0^  xr  /  A.  X* 

« 


1 

it  1 


(IV. 3. 5) 


and 


c  /  (<^+  L  )  k  xi  *  3L  +  JL.  }  (IV. 3. 6) 


while  the  decay  coefficient  X  has  the  form  given  by  Eq.  (11.76)  as 


X  -  (J+H  3} 

J,  -  a  A3k 

*  I  i  A 


(IV. 3. 7) 


iii)  Real  Gas  Analysis 

As  discussed  previously,  the  real  gas  behaviour  is  expressed  in  terms  of  two 
thermodynamic  properties,  namely  the  internal  energy  &  and  the  speed  of 
sound  factor  C  .  Thus,  before  attempting  to  solve  the  governing  equations, 
these  two  quantities  must  be  related  to  other  thermodynamic  variables  in  such 
a  manner  as  to  take  into  account  the  effects  of  compressibility  and  high  tem¬ 
perature.  To  this  end.  the  thermal  equation  of  state  may  be  expressed  in  the 
■form  of  a  perfect  gas  law  that  incorporates  a  compressibility  deviation  fac¬ 
tor,  2  ,  as  well  as  a  high  temperature  deviation  factor.  2r  .  Hence,  one 

v 
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may  write  an  equation  of  state  in  the  form 

9  =  zc  h  0 


(IV. 3.8) 


where 


0  *  JLE  (IV. 3. 9) 

with  R  being  the  gas  constant. 

One  may  then  think  of  determining  appropriate  forms  for  the  deviation  factors 

2C  and  £r  . 

The  compressibility  deviation  in  the  thermal  equation  of  a  perfect  gas  is  of¬ 
ten  taken  into  ucnsideration  by  means  of  semi-empirical  relations.  The  other 
properties,  like  internal  energy  and  specific  heats,  are  then  calculated  using 
the  concept  of  residual  properties  (Saad,  196?'*.  These  are  the  deviations  in 
the  real  gas  properties  from  those  of  the  perfect  gas  due  to  the  effect  of 
specific  volume  or  pressure  on  specific  heats.  For  example,  the  internal 
energy  may  be  expressed  as  €s  e\  er  where  e*  is  the  internal  energy  of  the 
perfect  gas  and  it  is  a  function  of  temperature  only,  cfr  is  the  residual  in¬ 
ternal  energy  and  it  is  a  function  of  both  temperature  and  specific  volume,  it 
is  given  by  Rozhdestrenskii  (1961)  as 

v 

=  /  [T(^=).-  P]  v  (IV. 3. 10) 

00  '  V 


where  V  is  the  specific  volume.  However,  the  thermal  equation  of  state  can 
be  put  in  the  form 


P  =  T  /  (V)  [l  -  ]  +  f  cv) 

1  L  Vr3J  * 


(IV. 3. 11) 
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where  the  functions  and  can  be  chosen  to  satisfy  equations  such  as  the 
Van  der  Waals'  equation  or  the  Beattie-Bridgeman  equation  and  C  is  a  constant. 
The  above  equation  is  found  to  fit  the  experimental  data  of  many  gases  to 
within  0.5%  accuracy  over  a  wide  range  of  pressures  and  temperatures. 

The  residual  internal  energy  may  now  be  calculated  by  differentiating  Eq.  (IV. 
3.11)  with  respect  to  temperature  and  then  carrying  out  the  integration  in  Eq. 
(IV. 3. 10).  This  yields 

V  C  r 

C  srl£  (  SL  Jv  _  f  f  Jv  (IV. 3. 12) 

r  ri  J  v  J  1 


The  residual  specific  heat  ,  difference  in  specific  heats,  Cp-Cv  and 

r  may  be  obtained  using  similar  relations  as  Eq.  (IV. 3. 10),  given  by  Saad 
(1969).  With  a  similar  procedure,  one  obtains  the  following  expressions: 


(IV. 3. 13) 


CP~Cv  = 


-T  (i+  v£)  ft 


T  ft  (l  -  )+  ft -i-  f, 


(IV. 3. 14) 


vr9 


and 


(IV. 3. 16) 


where  primes  indicate  differentiation  with  respect  to  V  . 

For  the  Beattie-Bridgeman  equation,  which  is  considered  here,  the  expressions 

for  f  and  f  and  their  derivatives  and  integrations  are 

*i  i 
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K-  £■  0-  Bo i  l 


*«•=£  i  (*  <*-*£> 


and 


/£  J*.=£[i*  $<*-&>]  ,/("-*<*-& 


(IV. 3. 16) 


where  A0  ,  <X  ,  0e  ,  b  and  C  are  constants  which  may  be  determined  experi¬ 
mentally  for  each  gas. 

To  define  the  caloric  equation  (  tf  *  ffr  )  completely,  an  expression  for 

* 

the  perfect  gas  internal  energy  e  should  be  given.  De  Broglie  suggested 
that  particles  on  the  atomic  scale  behave  in  a  wave-like  manner,  i.e.  associa¬ 
ted  with  each  particle  there  is  a  wave  with  a  certain  frequency.  The  Schroe- 
dinger  equation,  using  the  above  postulate,  describes  the  motion  of  the  parti¬ 
cle  in  terms  of  its  kinetic  energy  £  and  its  potential  energy  £  _  as 

Ki'rt  rot 

fOlloWSi 


2 

<7*  Y  +  g7r  2H  ((-.(:  )  ¥  = 

J*  K  Pot  '  T 


(IV. 3. 17) 


where  ^  and  m  are  the  displacement  and  the  mass  of  the  particle,  respective¬ 
ly,  while  ^  is  Planck's  constant,  and 


€  =  €  •+  € 

Kin  P„t 


Since  the  internal  energy  cf  the  perfect  gas  molecule  may  be  considered  to  be 
associated  with  translational,  rotational  and  vibrational  motion,  in  addition 
to  that  resulting  from  electronic  excitation,  the  Schroedinger  equation,  Eq. 
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(IV. 3. 17),  may  be  solved  for  each  mode  of  motion  independently  to  obtain  the 
various  available  energy  levels.  These  are  given  as 


and 


II 

b 

8** 

=  JL 

S7T* 

z  (f\  + 

(IV. 3. 18) 


where  n  ,  n,,  ,  n,  .  J  and  n  are  the  corresponding  quantum  numbers,  I  is 
the  moment  of  inertia  of  the  dumbell-shaped  molecule  and  -Q  is  the  frequency 
of  vibration  of  the  harmonic  oscillator;'  However,  for  a  system  of  W  partic¬ 
les,  each  single  molecule  can  exist  in  any  one  of  these  energy  levels,  i.e. 
there  is  a  certain  distribution  of  molecules  among  available  energy  levels. 

This  distribution  is  governed  by  the  thermodynamic  probability  factor  W  .  Ac¬ 
cording  to  the  Boltzmann  distribution  function  (Saad,  1969),  the  distribution 
of  particles  in  the  most  probable  case,  which  is  defined  by  W  =  0,  is  given 
by 


(IV. 3. 19) 


where  jj*  is  the  partition  function  and  is  defined  as 


Z  s 


(IV. 3. 20) 


gf;  is  the  degeneracy  (multiplicity)  of  energy  levels  with  the  same  value  of 
energy 

and  2  i _  with  4l  being  the  Boltzmann  constant. 

Tt 
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The  internal  energy  for  N  molecules  is  ?  s  aez-  N;  £(.  and  since  M(.  is  given 
by  Eq.  (IV, 3. 19),  the  caloric  equation  is  put  in  the  form 


a/4tz 

'  T»T  V 


(IV. 3. 21) 


which  is  evident  since 


~aZ  \  _  1 

~°t h  ~  hr* 


The  value  of  the  partition  function  may  be  evaluated  using  Eqs.  (IV. 3. 18)  and 
(IV. 3. 20),  after  the  summation  sign  is  approximated  by  an  integration  sign, 
since  the  separation  between  successive  energy  levels  is  small  (Lee  et  al., 
1963).  The  integration  yields: 


:hasij 


-2 


dr 


‘rot  = 


i  +  3  e  T  +  5  e 


-t - 


I  jr  (l  JL  0r  i_  ,drf 

l  K  *  3  T  +  15  {  T  ' 


+ - 


r<er 

)  T>6r 


and 


Vib 


-  3k 

e 


i  _  --wr 


:iv.3.22) 


where  Q  is  the  characteristic  temperature  of  rotation 

Y* 


Gy.  e 


-A 


zrllk 


-  152  - 


7^  is  the  characteristic  temperature  of  vibration 


T  s 
v  A 


Now,  according  to  the  type  of  the  molecules,  the  partition  function  and  the 
Internal  energy  of  the  perfect  gas  may  be  evaluated. 

For  a  monatomic  gas,  the  molecules  have  only  translational  motion,  accordingly 


v  (  z  ^ ^  T  )3/^ 

A* 

f  CIV. 3. 23) 


JL  A  NT 
z 


On  the  other  hand,  for  a  diatomic  gas,  the  molecules  have  translational,  rota¬ 
tional  and  vibrational  motion.  From  Eqs.  (IV. 3. 20)  and  (IV. 3. 21)  one  obtains 


V  /*7Tfn  )3/a  x  (—£. ) 

(  l3-  '  K  i-e-¥  ' 


(IV. 3. 24) 


It  is  to  be  noted  that  Eq.  (IV. 3. 24)  is  valid  for  T >  Qr  ,  which  is  the  case 
of  most  gases  at  T  >  100°K. 
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The  deviation  factor  2^  is  thus  obtained  from  Eq.  (IV. 3. 11)  and  Eqs.  (IV. 3. 16) 
in  the  following  form 


Zc  =  &  -  Bo  h  (4  -  [i  - 


^>3] 


-  Aq  A  ±(i_a*h) 


where 


e 


A*  *=■  -fix. 


Q  £Ta 


j  Bo  e  e0  4 


j  b'e  t  4 


(IV. 3. 25) 


and 

c  e  £.  & 

(QTJ 


The  high  temperature  deviation  factor,  2?r  ,  includes,  of  course,  the  effects 
of  both  dissociation  and  ionization  of  gases. 

During  the  process  of  dissociation,  the  composition  of  a  diatomic  gas  (\l  is 
given  by 


sz  (i  _  atj  )  A2  -f  1  otj  A 


(IV. 3. 26) 


where  aQ  is  the  degree  of  dissociation.  It  is  defined  as  the  ratio  of  the 
number  of  dissociated  it.,  lecules  to  the  initial  number  of  molecules  fi  .  The 
total  number  of  molecules  at  any  instant  is  given  from  Eq.  (IV. 3.26)  as 
(it  )  per  molecule  of  f\  .  Considering  each  constituent  as  a  perfect 
gas,  the  thermal  equation  of  the  dissociating  gas  may  be  written  in  the  form 
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P  =  (1  ^  )  S>  Rt  T 


(IV. 3. 27) 


where  is  the  gas  constant  of  f\z  . 

To  predict  the  equilibrium  composition  of  the  gas,  or  the  value  of  ag  as  a 
function  of  pressure  and  temperature,  the  law  of  mass  action,  Van’t  Hoff’s 
equation,  is  applied.  In  its  classical  form  it  is  given  as 


J  In  Kp 
AT  " 


R.T2 


(IV. 3. 20) 


where  Kp  is  the  equilibrium  constant,  and  <2p  is  the  heat  of  reaction  at  con¬ 
stant  pressure.  These  are  given  by 


Gp  =  ssz 


[&  ^t') p-  Li  V]  and  KP  *= 


7T 


V  f- 
r 


where  i)  is  the  number  of  molecules  and  is  the  enthalpy.  Subscripts  P 
and  r  denote  products  and  reactants,  respectively. 

The  partial  pressure  of  each  constituent  in  Eq.  (IV. 3. 26)  is  given  as 


1+  Xj 


(IV. 3. 29) 


_  i  -  <1 

i  +  *j 


Then  the  equilibrium  constant  may  be  calculated  from  Eqs.  (IV. 3. 283  and 
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(IV. 3.29)  as 


(IV. 3. 30) 


When  the  condition  of  statistical  equilibrium  is  applied  to  a  reacting  system, 
the  statistical  form  of  the  law  of  mass  action  is  obtained  in  the  form 


JL (  .ifgL  ) 
rtr  4 IT 


(IV. 3. 31) 


where  is  the  change  in  the  ground  state  energy  during  the  reaction. 

Equations  (IV. 3. 23)  and  (IV. 3. 24)  represent  the  partition  functions  of  the  re¬ 
actant  and  the  product  of  dissociation,  respectively .  Hence,  substituting  in 
Eq.  (IV. 3. 31),  one  obtains  the  equilibrium  constant  in  the  form  (Grossman, 
1969), 

0  3/a  3/2 

Or  )  T  [i-  exf(.  £)] 

(IV. 3. 32) 


where  is  the  dissociation  energy. 

The  caloric  equation  of  a  dissociating  gas  is  given  as 

^  —  (A-  -  +  2.  oL j  Gq  +  Gj  (IV.  3. 33) 

Differentiating  the  above  equation  with  respect  to  T  ,  one  obtains  the  speci¬ 
fic  heat  at  constant  volume,  Q ^  *  as 
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Cv  s  (4  -  ^  t-  (.2  (IV. 3. 34) 

The  value  of  (~“*L)V  is  obtained  by  differentiating  Eq.  (IV. 3. 30)  with  re¬ 
spect  to  temperature  and  using  Eq.  (IV. 3. 27)  to  substitute  for  P  .  After 
some  algebraic  manipulations,  one  obtains 


( 


-»T 


gk  /  *-<*  n  f~3  ,  ts_  72.  i 

T  -'La  T  "  t  (j?-Tv/r^  i) 


(IV. 3. 35) 


where 


rA 


*A 

■6 


Using  the  same  procedure,  the  specific  heat  at  constant  pressure,  Cp  ,  is  ob¬ 
tained  as 


:P  =  (i-Vj  )  q  t  c  ♦(*<;-  4,  ♦  ^)  Cg) 


with 


<-^)p - 


<sL  (±  _ 
zr 


V  (IV. 3. 36) 


Finally,  the  isentropic  compressibility,  or  the  speed  of  sound  factor,  f*  ,  is 
obtained  by  differentiating  Eq.  (IV. 3.27)  with  respect  to  ?  and  substituting 
in  Eq.  (IV. 3. 2).  It  is  found  to  equal 

r  -  c?  .  1 


(IV. 3. 37) 
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Thus,  Eq.  (IV. 3. 30)  and  (IV. 3. 31)  define  the  value  of  oCj  ,  and  Eqs.  (IV. 3. 34), 
(IV. 3. 36)  and  (IV. 3. 37)  determine  the  value  of  p  . 

Ionization  is  treated  very  similar  to  dissociation  from  the  point  of  view  of 
thermodynamics.  The  composition  of  the  ionized  gas  ft  at  any  temperature  is 
given  by 


A  =  A  +  oC(.  A+  +  ct;  If' 


(IV. 3. 36) 


where  a£{  is  the  degree  of  ionization  (subscript  L  =1,2  for  the  first  and 
second  ionizations,  respectively)  and  e~  denotes  a  free  electron.  It  is  de¬ 
fined  as  the  ratic  of  the  number  of  ionized  atoms  fp  to  the  original  number 
of  atoms  (\  .  The  thermal  equation  can  be  written  in  the  form 

P  =  (1  +  o( ■  )  f  HT  (IV. 3. 39) 


where  R  is  the  gas  constant  of  A  . 

The  equilibrium  constant  of  ionization  is  determined  by  substituting  from  Eq. 
(IV. 3. 23)  for  2  into  Eq.  (IV. 3. 31),  realizing  that  all  plasma  constituents, 
in  this  case,  can  be  treated  as  monatomic  gases.  Thus,  one  obtains 


<i 


TT  rn.e 

~T 


A 


3/2 

) 


(IV. 3. 40) 


where  €•  is  the  ionization 


td.  is  the  electronic 


U  is  the  electronic 

<e 

2  and  Zg  are  the 

A 

ft  and  A + 


energy 

mass 

partition  function  of  free  electrons  =  2 

partition  functions  of  electronic  excitation 
respectively . 


of 
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The  partial  pressure  of  each  constituent  is  evaluated  in  terms  of  <*■{  using 
Eq.  (IV. 3. 383,  as  previously  done  for  dissociation.  Substituting  in  Eq. 
(IV. 3. 283,  one  obtains 


a 


i  - 


(IV. 3. 413 


Equations  (IV. 3. 40)  and  (IV. 3. 413  give  the  equilibrium  composition  of  the 
ionized  gas  in  terms  of  its  pressure  and  temperature.  The  resulting  equation 
is  the  Saha  Equation  (Benson,  19673. 

The  internal  energy  of  the  ionized  gas  is  given  by 


and,  since  all  constituents  are  treated  as  monatomic  gases 


-  C‘ «■«,•)<+>  +  «f;  t.  ♦  etUc 


(IV. 3. 12) 


The  value  of  Cv  ,  Cp  and  P  are  obtained  in  a  similar  manner  as  in  dissocia¬ 
tion,  thus  one  gets 


and 


V  (IV. 3. 43) 


where 
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jo';  \ 

-»r  % 


r  (* 


while 


,2*i)  =  «L  (d-x*  )  (. 

K  -oT  '?  IT  L 


T;  s  J~i 


(IV. 3. 44) 


r  *  £l  , 


CV  1  +  Oi>  _  <¥• 


(IV. 3. 45) 


The  case  of  multiple  ionization  is  obtained  directly  by  modifying  the  set  of 
equations  of  single  ionization.  For  example,  the  Saha  Equation  may  be  exten¬ 
ded  (Gamble,  1963)  for  the  nth  electron  to  yield 


1-<X' 


=  t  A  (2rrmcA  a_ ,  Z-l ..  e*?(-  ) 

'  4*  ‘  .  P  '  Tt 


(IV. 46) 


Similarly,  Eqs.  (IV. 3. 42)  -  (IV. 3. 45)  may  also  be  extended  by  substituting 
to  replace  oc^  . 

Since  the  electron  motion  around  the  nucleus  is  very  complicated,  the  Schroe- 
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dinger  equation  is  not  solved  to  obtain  electronic  energy  levels.  However, 
since  the  separation  between  successive  energy  levels  is  high,  £  may  be 
expanded  in  the  form 


ze  -  %  *  %  e'& 


!  IV. 3. 47) 


where  both  the  degeneracies  and  the  energies  are  determined  experimen¬ 
tally  from  spectroscopic  data  (Rozhdestvenskii  et  al.,  1961). 

To  determine  the  electrqnic  internal  energy  and  specific  heat,  Eq.  (IV. 3. 47) 
is  differentiated  twice.  It  yields 


<  - 


- 


Twl  g  +.T 


(IV. 3. 46) 


-  few* 

( T*  _  £.)  Js! L  £  At 

v  r  4*  T3 


where  primes  denote  differentiation  with  respect  to  7*  . 

Substituting  from  Eq.  fIV.3.48)  in  Eq.  (IV. 3. 21),  the  internal  energy  is  ob¬ 
tained.  However,  by  differentiating  Eq.  (IV. 3. 21),  cv  is  obtained  in  terms 
of  2L  ,  Z*  and  Z*  ir  the  form 

*  C 


1  ( hiLzl£ 


RT  If 


(IV. 3. 49) 


The  main  difficulty  in  obtaining  solutions  with  the  series  expressed  in  Eq. 
(IV. 3. 47)  is  that  it  is  not  a  finite  sum.  The  terms  approach  infinity  as  the 
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energy  level  g,  approaches  a  constant  value,  namely  the  ionization  energy. 

However,  a  two-term  approximation  is  used  here  to  evaluate  the  thermodynamic 
properties. 

The  temperature  deviation  factor,  ,  is  then  obtained  from  Eq.  (IV. 3. 27) 
and  (IV. 3. 39)  for  a  double  ionized  gas  as 

2y  =  1  +  [if  Z  (1  +  1  <*a)]  (IV. 3. 50) 


where 


The  non-dimensional  equilibrium  constant  Kp  is  obtained  from  Eqs.  (IV. 3. 32) 

I  r 

and  (IV. 3. 46),  after  substituting  for  the  dimensionless  variables.  They  are 
given  as 


«5l.c  -^1 

4P  J  9  i/0,a 


[*-  c*r  (r 


&*)]  •*(-*£) 


and 


c,  j£ 


i.  5 


I'S 


e*P(-  &■) 


where 


2  A  #  (XjhaJl 


1-5 

) 


Jjg. 


and 
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The  pressure  j ^  is  in  atm.  and  the  temperature  T  is  in  °K. 

The  total  internal  energy  of  a  gas  may  be  expressed  as  the  summation  of  the 
contribution  of  various  effects,  i.e. 


.  * 

€  «=•  e  +  e 

C ©j*tp 


that  is  a  sum  of  the  ideal  gas  internal  energy  and  those  due  to  compressibi¬ 
lity,  dissociation,  ionization  and  electronic  excitation,  respectively. 

Therefore,  the  internal  energy  may  take  the  following  non-dimensional  form 


=  (i.  >  C«j-  +  ^  f  0j  H 

) 

*  C°1  ♦  <!-«,>  %  +  a 

+  «j  +  (1-^2) 

(IV. 3. 51) 


where  <T'  's  are  obtained  from  Eqs.  (IV. 3. 10),  (IV. 3. 23).  (IV. 3. 24),  (IV'. 3. 33) 
and  (IV. 3. 42)  as 


1B3  - 


6^  =  i.5  0 


s  2.5  9  +. 


<5-  =  -  *  C^A  ^)3  6?  (  i-e  Be*  k(o.5.  |Tk)} 

Covao  w  >■  ^ 


-  /»>  ^  <*-  #*^> 


where 


(IV. 3. 52) 


0  m  JL-  j  6  m  JL.  J  0{  e  -?- 

v  QTa  dGT*  QT« 


and  <s^f  =  0  . 

while  both  Be  and  are  given  in  Eqs.  (IV. 3. 47)  and  (IV. 3. 48). 

The  non-dimensional  speed  of  sound  factor  ,  which  is  defined  in  Eq.  (IV. 3. 2), 
can  be  put  in  the  form 


’sfiiP)  _  ^  s’*  (»  P  ) 

=TUfi  '  C 


(IV. 3. 53) 


The  value  of  #  is  given  by  Cp  /C^  .  which  are  obtained  from  Eqs.  (IV. 3. 34), 
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(IV. 3.36)  and  (IV. 3.43).  The  composition  of  the  gas  is  given  at  any  in¬ 
stant  by  Eqs.  (IV. 3.26)  and  (IV. 3. 38)  in  the  following  form 


♦  ^  +  «x  t  4^)  }]  *2«j 

*  4a  +  +  ** f«;„  +  4 y)]| 


(IV. 3. 55) 


are  given  in  Eq.  (IV. 3. 49).  while  ( 


are  expressed  in  Eqs.  CIV. 3. 35)  and  (IV. 3. 44). 

Since  various  energy  modes  are  excited  subsequently,  their  contribution  on 

(  HiflP )  may  be  expressed  as 

V  T 


_p{aps  rr  Q- 

T  Coi*tp  2  -y-  <3^-+- 


(IV. 3. 56) 


where  rfl  g  d  ,  1  and  2. 


The  term  )co**p  may  be  obtained  with  the  aid  of  Eqs.  (IV. 3. 9)  and 

(IV. 3.25)  which  yields 


,~z/?nP 


) 

Ccilf 


‘I  (i. 

9 


(IV. 3. 57) 


iv)  Solution 

a)  The  Self-Similar  Solution 

Since  self-similarity  reduces  appreciably  the  mathematical  complications  of 
the  solution,  it  has  always  been  the  first  step  in  seeking  a  solution  for 
blast  wave  problems.  The  utilized  equation  of  state  of  the  flowing  gas  should 
also  satisfy  the  conditions  of  self -similarity ,  i.e.  it  should  not  contain  any 
dimensional  constant  whose  dimensions  are  dependent  on  pressure.  As  previ¬ 
ously  suggested  by  Korobeinikov  et.al.  (1961)  and  Sedov  (195/),  its  general  form 
is 


e  =  p  <P(f) 


( IV  .,3.56) 
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where  <f>  is  a  function  of  density  only.  The  corresponding  thermal  equation 
may  be  derived  as  follows:  from  the  first  law  of  thermodynamics 


Tc/S  =  d  e  +  Pctv 


thus 


As  a  Jf+(2£.)  Ap]--±  JL  Af 

7*  jp  J  T  fz 

and  since  S  &  S  (Pj  *)  ,  then 

J  s  ~(^)  Jp+p»*\  A  f 

t 

Comparing  the  above  two  expressions  and  using  Eq.  (IV. 3. 57)  to  express  the 
differentials  of  the  internal  energy  in  terms  of  <p  and  some  thermodynamic 
identities  (Thompson,  1972),  one  obtains 


i  ( 


and 


W  (IV. 3. 59) 


where 


A  <P 
Jf 


The  above  two  equations  are  transformed  into  two  ordinary  differential  equa¬ 
tions  in  an  isentropic  process.  The  above  two  isentropic  relations  of  the 
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gas  obeying  Eq.  (IV. 3. 58)  are  obtained  by  integrating  them  as 

=C‘(S>  1 


and 


(IV. 3. 60) 


pep  e/pt-fJSL  )  =r  C2<S)  J 

Since  one  can  put  Cl  =  V'iCt  )  ,  the  general  thermal  equation  for  a  gas  obey 

ing  the  self-similarity  restrictions  is  obtained  as 


r  -  •*?  ( •*?({-  )) 


(IV. 3. 61) 


A  special,  but  important,  form  of  Eq.  (IV. 3.61)  may  be  obtained  if  is  ex¬ 
pressed  as  a  linear  function  of  its  arguments,  yielding 


r  =  b  p<p(f) 


(IV. 3. 62) 


where  Q  is  a  constant. 

The  Clausius  equation  of  state  is  suggested  to  represent  the  gas  behavior  at 
high  pressures  and  temperatures  (Thompson,  1972),  since  it  takes  into  account 
the  effect  of  the  volume  occupied  by  the  molecules.  As  it  satisfies  the  self- 
similar  conditions.,,  it  may  be  used  here.  This  equation  may  be  written  in  the 
form 


P  =  Rr./JL^  'j 
i-pf  J 


(IV. 3. S3) 
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where  p  is  a  constant. 

For  fluids  described  by  Eq.  (IV. 3. 63),  one  has  €  =  C(,T)  only,  just  as  the 
ideal  gas,  and  that  R  -  •  Thus 


e  =  _  R— . .  T 

if-  i 


(IV. 3. 64) 


with  Jf  being  a  mean  value  for  the  specific  heat  ratio. 

Comparing  Eqs.  (IV. 3.63)  and  (IV. 3. 64)  with  Eq.  (IV. 3. 58),  the  function  4>  is 
found  to  be 


<P 


__i _ L-f3*h 

k 


(IV. 3. 65) 


where 


P 


*r 


4 


Z3 


The  non-dimensional  speed  of  sound  factor,  C'  ,  is  found  by  substituting  in 
Eq.  (IV. 3.2)  from  Eq.  (IV. 3. 58)  and  Knowing  that 


then 


r »  jLf  p*'- p/j>t )  -  i-fW 


(IV. 3. 66) 


which  means  that  f  is  a  Function  of  density  only. 
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Substituting  from  Eq.  (IV. 3. 65)  into  Eq.  (IV. 3. 66),  one  gets 


r  = 


(IV. 3. 67) 


which  reduces  to  f  s  jC  for  a  perfect  gas  (  P  =0). 

The  governing  equations,  Eqs.  (IV. 3.1),  in  the  case  of  strong  explosion  condi¬ 
tion,  ^  =  0,  are  reduced  to 


<f-*)  J4.  +  h(4§  +  j£)  =o 


Jx  'Jx 

( f  X )  ^  4-  -r—  •  f 

(fx,7x*TJx  f 


=  o 


and 


(IV. 3. 60) 


The  decay  coefficient,  Eq.  (IV. 3. 7),  will  take  the  simple  form 


X  =  i*  i 


(IV. 3. 69) 


Equations  (IV. 3.68)  may  be  rewritten  in  their  autonomous  form  as 


dl 

dx 


■gam  +  .(ui)3 

k(f-xf  *9  — 


dk.,  -  _  Joz—  ( di  +j£) 
dot  (f-'X')  4*  X  ' 


(IV. 3. 70) 
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The  boundary  conditions  of  the  self-similar  flow  are  deduced  directly  from 
Eqs.  (IV. 3*3) .when  applying  the  condition  of  ^  =  0  and  consequently  <5^  =  0, 
one  obtains 


and 


ir  ( 


k„-l  f 

Kn  ' 


(IV. 3. 71) 


The  value  of  hn  is  evaluated  from  the  last  of  Eqs.  (IV. 3. 71)  when  is  sub¬ 
stituted  by  P  Q  cfo,  ,  then  after  some  algebraic  manipulations,  one  gets 

a  n 


A.rt  = - t±JL - 

S-d-f-,2 


and 


(IV. 3. 72) 


Xt  1 


) 


J 
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Since  the  boundary  conditions  depend  on  /3  ,  they  depend  on  the  initial  con¬ 
ditions  represented  by  and  hence,  the  whole  solution  is  expected  to  de¬ 
pend  on  the  initial  conditions.  However,  it  should  be  noted  here  that  the 
value  of  /3*  for  most  gases  (Saad,  1969)  is  in  the  order  of  0.002  at  the 
standard  atmospheric  conditions.  Thus,  its  effect  on  the  deviation  of  the 
gas,  and  consequently  on  the  whole  solution,  is  expected  to  be  negligible.  On 
the  other  hand,  for  moderate  values  of  initial  density  ,  the  problem  is 
greatly  affected  by  the  deviation  factor  since  it  Increases  with  . 

The  governing  equations  for  the  self-similar  solution,  Eqs.  (IV. 3. 70)  which 
are  subject  to  the  boundary  conditions,  Eqs.  (IV. 3. 72),  are  then  integrated 
numerically. 

» 

b)  The  Quasi-Self-Similar  Solution 

Once  the  self-similar  solution  is  established,  one  can  then  proceed  to  obtain 
the  non-self-similar  solution.  Before  attempting  to  solve  the  governing  equa¬ 
tions,  Eqs.  (IV. 3.1),  in  order  to  obta*''  the  flow  field  structure,  the  jump 
conditions  across  the  shock  have  to  be  determined.  The  Rankine-Hugoniot  equa¬ 
tion  is  found  by  substituting  for  o'  from  Eq.  (IV. 3.51)  into  the  last  of  Eqs. 
(IV. 3. 3)  as 


<^(Sn  -  e£ 


(IV. 3. 73) 


The  solution  procedure  may  be  outlined  as  follows 

a)  A  trial  value  of  hrt  is  assumed.  The  corresponding  value  of  the  per¬ 
fect  gas,  Eq.  (11.89),  at  the  same  front  Mach  number  may  be  used. 

b)  The  value  of  9n  is  calculated  from  Eq.  (IV. 3. 3)  as 


Jn 


) 


where 


Q 


may  be  taken  as 


since  the  undisturbed  medium  is  considered  as 
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a  perfect  gas. 

c)  Equation  (IV. 3. 8),  the  thermal  equation,  is  solved  numerically  with 

the  aid  of  Eqs.  (IV. 3. 25)  and  (IV. 3. 50),  for  the  temperature  .  The  false 
position  method  is  used  to  minimize  the  number  of  trials  needed  to  satisfy 
Eq.  (IV. 3. 8).  ' 

d)  The  values  of  and  Q ^  are  substituted  back  in  the  Rankine- 

Hugoniot  equation,  Eq.  (IV. 3. 73).  If  it  is  satisfied  tu  a  certain  permissible 
error  £  ,  the  assumed  value  of  is  taken  as  the  correct  one.  If  not,  an¬ 
other  trial  value  is  assumed  in  the  direction  which  reduces  the  error  and  the 
whole  precedure  is  repeated.  Here,  also,  the  false  position  method  may  be 
used  to  speed  up  the  convergence  of  the  solution. 

This  method  is  straight-forward  and  it  gives  for  any  values  of  0 ^  and  l~a  the 
properties  behind  the  shock  wave  up  to  any  desired  degree  of  accuracy.  In  ge¬ 
neral,  tse  equations  are  convergent  and  the  number  of  trials  do  not  exceed  5. 

Obtaining  the  boundary  conditions,  one  may  then  numerically  calculate  the  de¬ 
rivatives  of  the  gasdynamic  variables  with  respect  to  J  ,  which  are  equal  to 
their  values  at  the  front  as  stated  by  the  quasi-similar  technique.  The  New¬ 
ton  finite  difference  for  differentiation  is  U3ed  with  six  points  of  M  around 
the  value  of  M  at  which  the  field  is  calculated. 


The  computational  procedure  may  be  put  in  a  step-wise  manner  as  follows: 


a)  For  a  specified  shock  Mach  numbar  and  initial  condition,  the  boundary 
conditions  and  their  derivatives  are  calculated  numerically. 

b)  A  trial  value  of  A  is  assumed  and  the  integration  of  Eqs.  (IV. 3.1)  is 
carried  out  numerically  with  a  constant  step  size  U%  -  G.01.  The  temperature 

Q  is  calculated  at  every  step  in  the  integration  using  a  subroutine  which 
applies  the  false  position  method  to  solve  Eq.  (IV. 3, 8)  for  any  values  of  3 
and  h  . 


c)  The  integration  is  stopped  when  either  of  the 

is  encountered.  According  'o  the  lypu 

A 

decreased  or  increased. 


singularities  X  — or 
of  singularity,  is 


d)  When  the  value  of  is  sufficient  to  cairy  -.gt  tne  in .  ..-gi  ution  to  a 
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point  near  the  center  of  symmetry,  X  ==  ••  0.03,  the  velocity  profile  is  extra¬ 
polated  to  check  the  value  of  /  at  the  center,  f  ,  i.e.  .  The 

correct  solution  is  obtained  when  P(oj<4)  =  0,  Eq.  (11.112)  j  however,  if  fg<  € 
where  £  is  a  small  error,  the  corresponding  value  of  is  accepted.  If  f 

>  6  .  the  integration  is  repeated  for  values  of  \  around  this  value  and  a 
curve  between  ^  and  f  is  drawn,  the  correct  value  of  /!  is  that  correspon¬ 
ding  to  f  =  0  on  the  curve. 

0 

It  is  to  be  noted  that  the  quasi-similarity  approximation  does  not  conserve 
global  mass,  see  Sec.  (IV. 1),  i.e.  the  density  profiles  calculated  using  this 
method  do  not  satisfy  the  mass  integral,  Eq.  (IV. 3. 5).  On  the  other  hand,  the 
velocity  profiles  at  moderate  Mach  numbers  are  not  accurate  around  the  center 
of  symmetry  since  in  this  region  the  particle  velocity  should  be  negative  as 
the  particles  at  the  stage  of  expansion  return  back  to  their  original  position. 
However,  the  degree  of  accuracy  is  sufficient  for  a  good  qualitative  compari¬ 
son  between  solution  for  perfect  gases  and  real  gases  using  the  same  method. 

v)  Results  and  Conclusions 

Since  blast  waves  are  high  temperature  phenomena,  real  gas  effects  at  high  tem¬ 
peratures,  including  vibrational  excitation,  molecular  dissociation,  electronic 
excitation  and  ionization,  are  analyzed  to  determine  the  thermodynamic  proper¬ 
ties  of  gases  under  these  conditions.  The  condition  of  statistical  equilibrium 
between  various  species  in  a  chemical  reaction  is  used  to  evaluate  the  equili¬ 
brium  degree  of  excitation  of  the  reaction  processes  and  then,  the  thermodyna¬ 
mic  properties  are  calculated  considering  the  gas  to  be  a  mixture  of  perfect 
gases.  The  analysis  of  thermodynamic  properties  for  real  gases  is  suitable  for 
any  monatomic  or  diatomic  gas,  or  mixture  of  such  gases  as  in  the  case  of  air. 
However,  in  order  to  simplify  the  computations,  the  data  of  nitrogen  are  used 
in  this  application. 

The  effect  of  the  potential  energies  of  the  molecules,  which  may  be  divided  in¬ 
to  the  vibrational  energy,  energy  of  dissociation  and  energy  of  ionization,  on 
the  internal  energy  of  the  gas  are  represented  in  Fig.  IV. 3.1.  The  effect  of 
the  reacting  properties  on  the  deviation  factor  2^  ,  on  the  other  hand,  is  gi¬ 
ven  in  Fig.  IV. 3,2  as  a  function  of  pressure  and  temperature.  The  first  stage 


represents  the  increase  in  the  number  of  particles  due  to  dissociation,  thus 
indicating  the  degree  of  dissociation.  Similarly,  the  second  stage  shows 
twice  the  degree  of  ionization,  since  each  diatomic  molecule  splits  into  two 
mrnatomic  ones.  The  maximum  values  of  and  the  corresponding  minimum  va¬ 
lues  of  p  ,  Figs.  IV. 3. 3  and  IV. 3. 4,  occur  as  a  result  of  the  rapid  increase 
in  the  internal  energy  at  the  beginning  of  dissociation  and  every  successive 
ionization. 


From  the  physical  point  of  view,  the  function  of  density,  <p  ,  in  the  self- 
si.nj.lar  solution  should  be  selected  to  satisfy  one  of  the  well-known  equa¬ 
tions  of  state  suitable  to  be  applied  in  the  range  of  conditions  of  self¬ 
similar  flow.  Thus,  the  Clausius  equation  of  state  is  used  to  represent  the 


non- perfect  gas  deviation  due  to  compressibility  at  high  pressures,  in  terms 
of  ,9*r  which  depends  on  A  .  However,  since  both  <p  and  the  boundary  con¬ 
ditions  depend  on  fa  ,  contrary  to  the  casB  of  perfect  gases,  the  whole  so¬ 
lution  should  be  carried  out  for  specified  values  of  Pa  .  The  results  re¬ 
presented  in  Figs.  IV. 3. 5  -  IV. 3. 8  are  obtained  for  values  of  -  O.Q, 

0.016  and  0.032.  These  values  are  corresponding  to  high  values  of  fa  at 
which  the  compressibility  effect  due  to  high  pressures  overcomes  the  effect 
of  high  temperatures  that  reduces  such  effect.  Of  course,  at  these  high  tem¬ 
peratures,  the  effect  of  compressibility  is  neglected  at  normal  densities. 

The  results  are  calculated  for  J  =  2  and  if  =  1.4. 

Oue  to  the  complicated  form  of  the  equation  of  state  in  t he  case  of  rea' 
ses,  the  variati--*  of  the  gasdynamlc  properties  acrocs  the  shock  wave,  wh.ch 
represent  the  boundary  conditions  of  the  problem,  are  evaluated  numerically 
by  solving  the  Rankine-Hugoi.iot  equ>.  ion  with  the  caloric  equation  of  state 
using  the  false  position  method  to  solve  a  system  of  algebraic  equations. 

The  density  ratio  increases  with  the  excitation  of  successive  energy  mu- 

t 

des,  Fig.  IV.3.SJ.  From  cq.  (IV.  3. 3}  for  strong  shoe*,  waves,  sf  s-  jL(i_  JL )  . 

n  3  Kn 

Thu-..,  db  e-'  increases  due  to  the  effect  of  the  potential  energy,  kn  also 
increases.  Gn  the  other'  hand.,  the  perfect  gas  density  ratio  for  strong  shock 

waves,  which  is  given  by  kn  .  may  help  in  explaining  the  sequence 

If-  1 

of  events  of  Fig.  iV.3.9.  Since  <jf  Increases  due  to  viorational  excitation 
and  reaches  a  minimum  at  the  beginning  of  dissociatic  ,,  hn  increases  to  reach 
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a  maximum  at  the  same  point.  Also,  the  decrease  of  hn  with  Mach  number  during 
ionization  is  related  to  the  increase  in  %  ,  since  at  this  region  Cv  increa¬ 
ses  due  to  the  increase  in  the  number  of  particles  which  overcomes  the  effect 
of  the  potential  energy  (Zel’dovich  et  al.,  1966). 

The  temperature  ratios  decrease  with  respect  to  those  of  perfect  gases  due  to 
the  expenditure  of  energy  as  potential  energy,  as  well  as  the  increase  in  the 
number  of  particles.  The  high  values  of  potential  energy  as  compared  to  the 
translational  energy  and  the  increase  in  decrease  the  rate  at  which  the 
temperature  ratio  increases  with  Mach  number  as  detected  from  Fig.  IV. 3. 10. 

Figure  IV. 3. 11  shows  that  the  pressure  ratio  is  not  affected  appreciably  by 
the  real  gas  properties.  This  can  be  qualitatively  seen  from  Eqs.  (IV. 3. 3) 
which  indicate  that  the  pressure  ratio  varies  as  »  which  is  not  affected 

by  the  real  gas  properties  inside  the  blast  wave.  However,  since  |\n  increases 
for  real  gases,  the  pressure  ratio  also  increases. 

Actually,  the  density  ratio  curve  has  a  number  of  maxima  which  occur  at  the  be¬ 
ginning  uf  the  successive  ionizations  like  that  occurring  at  the  beginning  of 
dissociation.  After  the  occurrence  of  total  ionization  the  gas  transforms  into 
a  mixture  of  monatomic  particles  which  have  a  constant  value  of  If  =  1.667. 

Thus,  asymptotically  approaches  the  value  of  four. 

The  effect  of  initial  pressure  on  density  ratio  and  temperature  ratio  across 
the  shock  wave  is  shown  in  Figs.  IV. 3. 12  and  IV. 3. 13,  respectively.  At  a  cer¬ 
tain  Mach  number,  as  the  ambient  pressure  increases,  the  temperature  ratio  in¬ 
creases  while  the  density  ratio  decreases. 

The  flow  fields  generated  by  blast  waves  at  various  values  of  front  Mach  numbers 
are  evaluated  by  integrating  the  system  of  equations,  Eqs.  (IV. 3.1),  with  f 
expressed  by  Eq.  (IV. 3. 53).  Using  the  quasi-similarity  approximation,  the  deri¬ 
vatives,  with  respect  to  J  ,  are  substituted  for  by  their  corresponding  values 
at  the  front.  The  integration  is  carried  out  numerically  using  the  Runge-Kutta 
fourth  order  method  with  a  step  size  -  0.01.  The  computations  are  performed 
for  the  spherical  wave  only,  j  =  2. 

The  decay  coefficient  %  >  which  is  calculated  using  an  iterative  procedure  to 
satisfy  the  compatcbility  condition  at  the  center  of  synmetry,  decreases  with 


the  excitation  of  the  successiv/e  energy  modes.  Fig.  IV. 3. 14.  Although  this  de¬ 
crease  seems  to  be  slight,  its  effect  is  rather  important  since  the  field  is 
very  sensitive  to  the  value  of  X  ard  any  change  in  its  third  decimal  value  af¬ 
fects  the  field  appreciably.  This  weak  dependence  of  the  decay  coefficient  on 
the  field  was  previously  observed  by  Kamel  (1973b).  However,  the  slight  de¬ 
crease  in  value  of  X  roay  be  due  to  the  expenditure  of  energy  in  the  excitation 
processes  which  can  be  represented  as  an  energy  sink.  1  he  deviation  in  X  from 
its  corresponding  value  of  the  perfect  gas  decreases,  approaching  zero  at 
fi  =  2.5  where  temperature  inside  the  field  is  not  sufficiently  high  to  excite 
any  additional  energy  modes. 

It  is  logical  that  temperature  distributions  inside  the  field  are  greatly  affec¬ 
ted  by  real  gas  behaviour  since  a  great  amount  of  the  blast  energy  is  used  as  a 
potential  energy,  beside  the  translational  energy.  It  is  to  be  noted  that  when 
comparing  the  temperature  distribution  of  the  perfect  gas  field  with  that  of 
the  real  gas,  other  effects  beside  the  energy  of  excitation  should  be  taken  into 
consideration.  For  example,  when  the  vibrationa)  motion  of  the  molecules  are 
excited,  an  amount  of  energy  is  expended  which  would  decrease  the  temperature. 
But,  on  the  other  hand,  the  excitation  of  vibrational  motion  of  the  atoms  rela¬ 
tive  to  each  other  increases  the  volume  occupied  by  the  molecules  and,  conse¬ 
quently,  increases  the  temperature.  The  predominant  factor  of  these  two  effects 
depends  on  the  degree  of  excitation  of  the  vibrational  motion.  The  resultant 
effect  may  be  observed  in  the  field  of  hi  =  3  and  5  in  Figs.  IV. 3. 15  and  IV. 3. 16. 
Figure  IV. 3. 17  shows  a  comparisor  between  the  temperature  profiles  for  the  per¬ 
fect  gas,  with  different  ,  and  for  the  real  gas,  while  Fig.  IV. 3. 18  shows  the 
effect  of  initial  pressure  on  the  temperature  profiles. 

However,  in  the  process  of  dissociation  and  ionization,  the  heat  of  reaction  is 
very  high  and,  thus,  the  decrease  in  temperature  is  always  appreciable  and  the 
portions  of  lower  gradients  in  the  temperature  curves.  Figs.  IV. 3. 15  through 
IV. 3. 18  represent  the  successive  reacting  processes. 

The  density  profiles.  Figs.  IV. 3. 19  through  j.V.3.22,  are  related  directly  to  the 
temperature  profiles.  These  figures  indicate  that  the  density  increases  behind 
the  shcc.N  front  and  in  the  vicinity  of  the  center  of  symmetry. 
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In  contrast  to  the  large  deviations  in  temperature  and  density  profiles  from 
those  for  perfect  gases,  the  pressure  and  particle  velocity  profiles  are  found 
to  be  less  sensitive,  Figs.  IV. 3.23  through  IV. 3. 25.  This  is  a  direct  conse¬ 
quence  of  the  fact  that  both  pressure  and  particle  velocity  are  related  to  the 
particle  translational  energy.  Although  the  direct  proportionality  between 
the  particle  translational  energy  and  the  temperature  would  suggest  that  the 
former  should  decrease  sharply,  according  to  the  behaviour  of  the  latter  given 
by  Figs.  IV. 3. 15  through  IV. 3. 18,  the  large  increase  in  i.ie  number  of  partic¬ 
les  due  to  dissociation  and  ionization,  tend  to  offset  this  effect,  hence  the 
insensitivity  depicted  in  Figs.  IV. 3. 23  through  IV. 3. 25. 

Finally,  the  changes  in  f  and  C ^  within  the  flow  field  are  illustrated  in 
Figs.  IV. 3. 26  through  IV. 3.29,  where  their  behaviour,  as  in  the  case  of  tempe¬ 
rature,  manifests  the  different  modes  of  excitation  within  the  blast  wave  struc¬ 
ture. 
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Fig.  IV. 3.1  Real  gas  effect  on  internal  energy  at  different  values  of  the 
ambient  pressure  for  4  ~  2. 

Fig.  IV. 3. 2  Real  gas  deviation  factor  F  as  a  function  of  temperature  at 
different  values  of  the  ambient  pressure  for  4  =2. 

Fig.  IV. 3. 3  Effect  of  vibrational  exitation,  dissociation  and  ionization  on 
P  at  different  values  of  the  ambient  pressure  for  j  =2. 

Fig.  IV. 3. 4  Effect  of  dissociation  and  ionization  on  at  different  va¬ 
lues  of  the  ambient  pressure  for  i  =2. 

Fig.  IV. 3. 5  Self-similar  density  profiles  at  different  values  of  the  para¬ 
meter  p*  with  S'  =  1.4  while  d  =2. 

t 

Fig.  IV. 3. 6  Self-similar  pressure  profiles  at  different  values  of  the  para¬ 
meter  p*  with  =1.4  while  4  =2. 

Fig.  IV. 3. 7  Self-similar  velocity  profiles  at  different  values  of  the  para¬ 
meter  p*  with  8”  =  1.4  while  4  =  2. 

Fig.  IV. 3. 8  Self-similar  temperature  profiles  at  different  values  of  the 

parameter  p  with  if  =  1.4  while  J  =2. 

Fig.  IV. 3. 9  Real  gas  effects  on  density  ratio  across  a  shock  wave. 

Fig,  IV. 3. 10  Real  gas  effects  on  temperature  ratio  across  a  shock  wave. 

Fig.  IV. 3. 11  Comparison  between  real  and  perfect  gas  pressure  ratios  across 
a  shock  wave. 

Fig,  IV. 3. 12  Effect  of  initial  pressure  on  the  density  ratio  across  a  shock 


wave . 


Fig.  IV. 3. 13  Effect  of  initial  pressure  on  the  temperature  ratio  across  a 
shock  wave. 

i 

Fig.  IV. 3. 14  Real  and  perfect  gas  decay  coefficients  as  a  function  of  Jj  for 

d  =  2. 

Fig.  IV. 3. 15  Temperature  distribution  for  various  excitation  of  internal 
energy  modes  at  H  =  3  for  j  =2. 

Fig.  IV. 3. 16  Temperature  distribution  for  various  excitation  of  internal 
energy  modes  at  H  =  5  for  d  =2. 

Fig.  IV. 3. 17  Comparison  between  the  temperature  profiles  for  the  perfect  gas, 
with  different  values  of  ,  and  for  the  real  gas  at  H  c  7  for 

i  =  2. 

Fig.  IV. 3. 18  Initial  pressure  effect  on  temperature  distribution  at  M  ■*  IQ 

for  j  =>  ,2. 

» 

tig.  IV. 3. 19  Real  gas  effect  on  density  profiles  at  fl  =  3  for  d  =2. 

Fig.  IV. 3.20  Density  profiles  for  various  excitation  modes  at  *1  =  5  for 

j  =  2. 

Fig.  IV. 3. 21  Comparison  between  the  density  profiles  for  the  perfect  gas, 
with  different  values  of  X  ,  and  for  the  real  gas  at  H  -  7 
for  j  =2. 

fig.  IV. 3. 22  Initial  pressure  effect  on  dens? r,  profiles  at  M  =  10  for 

i  =2. 

Fig.  IV. 3. 23  Velocity  and  pressure  distributions  for  both  perfect  and  real 
gases  in  a  spherical  blast  wave  at  Kl  =5. 
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Fig.  IV. 3.24  Real  gas  effect  on  velocity  and  pressure  distributions  as  com¬ 
pared  with  those  of  perfect  gas  with  different  values  of  at 
H  =  7  for  a  spherical  wave. 

I 

Fig.  IV. 3.25  Comparison  between  the  velocity  and  pressure  profiles  for  the 
perfect  gas  and  for  the  real  gas,  with  different  initial  pres¬ 
sures,  at  rj  =  10  for  i  =  2. 


Fig.  IV. 3.26 


Distribution  of  p  for  various  excitation  modes  at  M  =  3  for 

j  =  2. 


Fig.  IV. 3. 27  Distribution  of  P  for  various  excitation  modes  at  fi  =  5  for 
d  =  2. 


Fig.  IV. 3.28  Distribution  of  C '  with  successive  excitation  at  M  =  3  for 

j  =  2. 

Fig.  IV. 3.29  Distribution  of  with  successive  excitation  at  H  =  5  for 

J  -  2. 
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XV. 4.  BLAST  WAVES  IH  A  DETONATING  MEDIUM* 


i)  Introduction 

In  this  section,  a  study  oP  the  gasdynamic  aspects  accompanying  the  propagation 
of  e  detonation  wave  in  a  uniform  combustible  medium  is  presented. 

According  to  the  Chapman-Jouguet  theory  (Cppenhe:'.m,  1968},  a  detonation  wave, 
once  formed,  will  propagate  at  a  constant  velocity  independent  of  the  details 
at  its  initiation  processes,  and  dependent  only  cn  the  properties  of  the  explo¬ 
sive  medium  in  ths  front  of  tha  wave. 

The  simplest  model  to  use  as  a  basis  for  an  initiation  theory  would  be  the  in- 

i 

stantaneous  formation  of  a  reactive  blast  wave  front. 

Since  methane  is  the  major  constituent  of  natural  gas,  the  most  predominant ly 
used  explosive  gas  mixture,  the  detonating  medium  in  this  application  is  as¬ 
sumed  to  be  a  stoichiometric  methane-air  mixture.  The  detonation  medium  itself 
is  taken  as  a  perfect  gas  with  a  mean  specific  naat  ratio,  ,  that  is  capable 
of  producing  an  amount  of  exothermic  energy,  ,  per  unit  mass  of  the  mixture. 

•The  solution  obtained  here  is  for  the  non-sslf-similar  flow  fields  that  result 
from  instantaneous  deposition  of  energy  at  a  point,  along  a  line  or  at  a  plane 
in  an  inviscid  medium.  It  traces  the  decay  of  the  datonation  front  from  its 
initial  sslf-similar  limit,  corresponding  to  the  adiabatic  point  explosion  to 
its  Chapman-Jouguat  condition. 

■i i 5  Problem  Formulation 

Tha  basic  equations  governing  the  p-ebiera  can  be  specified  directly  from  Eqs. 

( 1 1 . 2  C  jf ,  (11.22}  and  (11.24)  by  equating  all  *s  as  well  as  fa  's  to  zero, 
since  wg  ‘-ave  a  sourceless  flow  field  and  a  uniform  ambient  medium.  These 
equations  are 


"■This  application  is  based  c-n  Ahdei-Rauuf 


11962} . 


-  211 


™  $■+<*■*>&  +  ■« 
+  *£■  ;*  c0 


n 


(IV. 4.1) 


a  s  ai  _  AJ 

-»s 


<;-«  *|  +  <s<5|  +  j£)«« 


j 


The  boundary  conditions  of  the  problem  are  given  by  the  Hugoniot  relations, 
Eqs.  (11.93),  (11.54)  and  (11.85),  which  may  be  rewritten  as 


and 


(IV. 4. 2) 

o 
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The  shGcK  strength  at  the  Cnapman-Jouguet  condition  of  overdriven  detonation 
is  given  by  Eq.  (11.94)  as 

*CU  =  i  +  *]  1)5+  if-  i  (IV.  4,3! 

The  mass  and  energy  integrals  for  the  present  problem,  Eqs.  (II. 64 i  arid  (Il.bUi 
respectively,  are  given  by 

=  (  kx*  S  (T 

a  J+  i 
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(IV. 4. 5) 


The  mass  integral,  Eq.  (IV. 4. 4),  is  used  here  as  a  check  on  the  density  pro¬ 
files,  while  the  energy  integral,  Eq.  (IV. 4. 55,  is  used  for  determining  the 

non-dimensional  radius,  J  . 

Inc'  decay  coefficient,  A  ,  takes  its  form  given  by  Eq.  (11.97) 


(i+1) 


v<r-  i) 


0, 


-y 


d  J| 
d  * 


(IV. 4. 6) 


which  for  the  strong  self-similar  blast  wavs  (  t/  *>  0),  reduces  to 


'A  =(Ui) 


i IV. 4.7) 


Finally,  trie  compatability  condition  which  putt  be  satisfied  at  the  center  of 
symmetry  by  the  correct  solution  is  given  by  Eq.  (TI.H2),  namely 


/  (o>y)  =  o 


(iv.  4., a ) 


applying  the  concept  of  the  q.„asi  similar  technique  cn  the  governing  equal  i...ns. 
'viilS  tu  tnc  following  uruiud-y  differential  equations,  giueu  picwic 


1  K&j 


i  .  -  l  .  , 

v’ 
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by  Eqa.  (Ij.1,23)  -  (111,25),  when  all  's  vanish 

(/-  X )  dk.  +  k  (dL  +  ]£.)+  %k  ft  zs  o 

6X  dX  PC 

(f-z)  M.  +  ±  dl  +  6  so 

dX  K  dX 


and 


<*-*>£+'  *(&+ 


/v 


ss  O 


where 


ft  m  — 


0  m 


o.  S 


and 


Ca  e  %  ^ 

Xf»  +  n  ~ 


with 


Equations  (IV. 4. 9)  thus  yield 
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B 


and 


21  = 

Ax  k  (f~x)*  -  *3 

dk  =  jJl  (dL+  3£  +  %A) 

4%  (p~x)K  7x  x 


(IV. 4. 10) 


JL  «  _  h  [(f-x)J£  +  A /’s'] 


P 


As  discussed  previously  in  Chapter  III,  there  are  two  different  singularities, 
in  addition  to  the  singularity  at  the  center  of  explosion,  which  can  be  direct¬ 
ly  observed  from  Eqs.  (IV. 4. 10)  These  singularities  are  located  at  f  •  X>  and 
K  (f  -  X)  31  IfS  ,  In  terms  of  the  reduced  variables,  f  and  Z  .  these  sin¬ 
gularities  are  given  by 

f"  m  1-0  and  »  ("1-F")  (IV. 4. 11) 


iii)  Solution 

Equations  (IV. 4. 10)  can  be  integrated  numerically  using  Runge-Kutta  fourth  order 
not  hod ,  subject  to  the  boundary  conditions  given  by  Eqs.  (IV. 4. 2).  An  iterative 
procedure  must  be  used  to  find  the  correct  value  of  %  which  satisfies  Eq. 

IIV. 4. 8)...  Since  we  cannot  proceed  with  the  numerical  integration  till  the  cen¬ 
ter  of  syrrunstry,  due  to  the  existence  of  the  singularity  at  the  center,  the  so¬ 
lution  is  based  on  matching  the  numerical  solution,  from  X  -  1  to  %  ®  X  , 
with  an  asymptotic  formula  fur  the  non-dimensional  vc)nr,ioy,  f  .  However,  it 
can  be  assumed  that  the  velocity  distribution  near  the  center  of  symmetry  is 
linear,  provided  that  the  value  of  X  is  very  small.  The  value  of  X  nmust  be 
found  a  priori  for  each  geometrical  symmetry. 
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The  (4  domain  in  which  we  are  Interested  (  9  <  !j  <  )  la  divided  Into 

small  steps,  each  step  ( 4  if  ]  equals  0.005.  The  choice  of  the  values  of 
and  4^  Is  based  on  the  grid-independent  solution. 

If  the  solution  at  ^  Is  Known,  then  the  solution  at  54  ■ 

(  a  %  +  4*4  )  is  obtained  as  follows! 

I 

a)  At  a  certain  value  of  •  54, vi  ,  the  boundary  conditions  fn  ,  h„ 

and  <3n  (at  %  ■  1)  can  be  calculated  from  Eqs.  (IV. 4. 2). 

b)  An  assumed  value  of  the  decoy  coefficient  is  then  used  to  inte¬ 

grate  Eqs.  (IV. 4. 10)  numerically.  The  integration  is  stopped  at  %  ■  #  . 

c)  Applying  the  linear  relation  between  the  non-dimensional  velocity  f 
and  X  from  9C  •  %  to  X  •  0,  the  velocity  at  the  center  of  symmetry  can  be 
obtained  by  extrapolation. 

d)  If  the  absolute  value  for  the  velocity  at  the  center  of  symmetry,  which 

should  identically  be  equal  to  zero,  is  greatsr  than  ths  required  accuracy, 
one  may  repeat  the  above  procedure  by  changing  the  assumed  by  lX  until 

the  required  accuracy  la  obtained. 

e)  With  the  correct  value  of  determined,  the  flow  field  parameters 

will  be  obtained,  and  then  the  values  of  the  mass  integral,  ,  and  the  ener¬ 
gy  integral,  0^  ,  can  be  calculated  using  Romberg  procedure. 

i  .  i 

f)  The  whole  procedure  is  repeated  for  other  values  of  54  until  the  solu¬ 
tion  is  obtained  for  the  whole  range  of  54  . 

It  should  be  noted  that  forcing  the  solution  to  proceed  to  too  small  a  value 
of  X  will  throw  the  solution  into  the  singularity  X  *  0.  Also,  during  the 
trials  of  obtaining  the  correct  X  »  caution  should  be  exercised  lest  the  solu¬ 
tion  throw  into  one  or  the  other  of  the  singularities  given  by  Eq.  (IV. 4. 11). 


IV.  Results  and  Conclusions 

In  order  to  obtain  a  numerical  solution,  it  has  been  assumed  that  the  detona¬ 
ting  medium  is  a  stoichiometric  methane-air  mixture.  Therefore,  the  value  of 
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$  for  «uch  a  mixture  la  20.0231  and  S  *  1.3  (Kamel  at  al.,  1979).  For  the 
specific  case  presented  here,  Eq.  (IV. 4. 3)  yields  -  0.033786. 

Figures  IV. 4. la,  b  and  c  show  the  pressure  profiles  for  planar,  cylindrical 
and  spherical  geometries  with  S’  ■  1.3.  At  any  specified  value  of  y  ,  the 
pressure  immediately  behind  the  wave  front  is  a  maximum  and  falls  off  to  a 
nearly  constant  value  near  the  center  of  symmetry.  The  fall-off  in  pressure 
is  greatest  for  spherical  waves,  with  its  greatest  value  of  the  decay  coeffi¬ 
cient  'X  ,  and  least  for  the  planar  waves.  As  in  the  case  of  adiabatic  point 
•  explosions  (Korobeinikov  et  al.,  1969),  the  pressure  ratio  at  the  center  de- 

i 

creases  first  to  a  minimum  before  ascending  to  its  Chapman -Jouguet  detonation 
value. 

Figures  IV. 4. 2a,  b  and  c  and  IV. 4. 3a,  b  and  c  show  the  density  and  temperature 
profiles,  respectively.  These  sets  of  curves  are  bound  to  the  right  by  the 
self -similar  solution  for  the  adiabatic  point  explosion.  It  is  evident  that 
the  density  at  the  center  is  always  zero,  corresponding  to  infinite  tempera¬ 
ture,  a  characteristic  property  of  the  solution  obtained  for  an  invisoid,  non¬ 
conducting  and  non-radiating  gas.  As  shown  in  Figs.  IV. 4. 2,  nearly  all  of  the 
mass  of  the  gas  engulfed  by  the  blast  is  concentrated  close  to  the  wave  front 
for  the  self-similar  case  and  is  then  distributed  gradually  in  the  flow  field 
for  >  0,  until  the  Chapman-Jouguet  state  where  the  distribution  of  the 
mass  is  more  uniform.  It  is  seen  also  that  the  fall-off  in  density  1b  greatest 
for  spherical  waves,  the  same  tendency  in  pressure  profiles,  and  least  for 
planar  waves. 

The  particle  velocites, 'as  shown  in  Figs.  IV. 4. 4a,  b  and  c.  have  a  maximum  va¬ 
lue  immediately  behind  the  wave  front  and  decrease  to  zero  at  the  center  of 
symmetry.  At  any  given  %  and  Jj  ,  the  spehrical  flow  velocity  is  the  lowest. 
With  larger  values  of  X  ,  and  for  planar  waves  the  velocity  increases. 

The  integral  curves  for  different  values  of  y  ranging  from  ■  0  to  *d  ■  0.033, 
are  projected  on  the  f  -  Z  phase  plane,  for  4  ■  0,  1  and  2  in  Figs.  IV. 4. 5a, 
b  and  c,  respectively.  Integral  curves  for  if  ■  0  represent  the  solution  at 
the  initial  instance,  whils  the  initiation  energy,  Ej  ,  is  still  the  predomi¬ 
nant  parameter  governing  the  flow  field.  They  are,  of  course,  identical  to 


those  for  the  self-similar  point  explosions  in  an  inert  gas  (Sedov,  1957).  As 
the  front  decays  and,  concomitantly,  the  role  of  chemical  energy  becomes  more 
prominent,  the  value  of  $  increases  until  it,  finally,  asymptotically  reaches 
the  steady  Chapman- Jouguet  value  characterizing  the  explosion.  Each  integral 
curve  starts  from  a  particular  point  on  the  Hugoniot  curve,  specified  by  the 
conditions  corresponding  to  %  -  i  .  It  can  be  directly  observed  from  Figs. 

IV. 4. 5  that  the  solutions  are  free  from  singularities  since  the  locii  of  the 
singularity  z  and  the  singularity  F  *  1  do  not  intersect  with  any 

of  the  integral  curves.  Since  a  singularity  does  exist  at  the  C-J  condition, 
the  corresponding  integral  curve  is  not  included.  One,  however,  may  approach 
that  condition  sufficiently  to  get  a  good  indication  of  the  quantitative  and 
qualitative  nature  of  the  solution. 

The  mass  integral,  (3^  ,  and  the  energy  integral,  ,  are  depicted  in  Figs. 
IV. 4. 6  and  IV. 4. 7,  respectively,  for  different  values  of  j  .  These  curves 
start  from  thQir  respective  values  for  the  inert  self-similar  point  explosions 
at  y  °  0  ;.nd  end  at  i/  «  0.033  which  nearly  equals  .  The  mass  integral 

is  used  to  check  the  accuracy  of  ths  solution  and  its  exact  value  is  given  by 
Eq.  IV. 4. 6  that  some  error  is  obtained  in  the  mass  integral  which  reflects  in¬ 
accuracy  in  the  density  distribution  within  the  flow  field.  This  error,  how¬ 
ever,  is  expected  due  to  the  quasi-similar  approximation  which  causes  the  solu¬ 
tion  to  be  accurate  near  the  self -similar  range  and  just  behind  the  wave  front. 
The  energy  integral  may  be  used  to  iterate  on  the  decay  coefficient,  which  is 
another  method  to  obtain  the  correct  %  ,  but  it  is  used  here  to  obtain  the 
non-dimensional  wave  radius,  J  ,  using  Eq.  IV. 4. 5. 

i he  decay  coefficients,  %  ,  are  shown  as  functions  of  tj  in  Fig.  IV. 4. 8.  An 
iterative  procedure  was  performed  to  obtain  the  correct  value  of  A  which  sa¬ 
tisfies  the  condition  of  zero  particle  velocity  at  the  center  of  symmetry.  Any 
assumed  value  of  *X  less  or  greater  than  the  correct  one  results  in  a  numerical 
instability .  This  occurs  due  to  the  two  different  singular  points  in  the  flow 
field,  given  by  Eq.  (IV. 4. 11),  which  would  pull  the  integral  curve  to  either  of 
them.  It  was  observed  that  the  numerical  instability,  due  to  these  singulari¬ 
ties,  may  occur  in  two  different  forms.  The  dependent  variables  p  and  ft  may 
take  negative  values  or  may  take  values  greater  than  the  preceeding  ones  during 


the  numerical  integration  of  the  governing  equations.  Therefore,  it  is  recom¬ 
mended  to  use  the  conditional  "IF"  statement  inside  the  integration  procedure 
in  the  computer  program  to  avoid  the  singular  points  and  save  computational 
time  during  seeking  the  correct  \  .  It  was  also  found  that  a  deviation  in 
the  order  of  1  x  10  *  in  the  correct  value  of  may  lead  to  the  domain  of 

singularities.  Because  of  the  singularity  F  =1.0,  the  numerical  integration 

* 

of  the  governing  equations  was  stopped  at  X  =  0.02  for  planar  and  cylindrical 
waves  and  at  %  •  0.04  for  spherical  waves.  A  linear  extrapolation  was  then 
applied  to  f  in  order  to  obtain  f  (0,  !j  ).  The  solution  was  considered  cor¬ 
rect  when  i  P(©jtO|  x  10  .  In  order  to  achieve  this  accuracy,  the  improve- 

A  “4  "6 

ment  in  the  value  of  a  was  in  the  order  of  1  x  10  for  planar  waves,  1  x  10 

.0 

for  cylindrical  waves  and  1  x  10  for  spherical  ones. 

The  non-dimensional  radius,  J  ,  is  plotted  as  a  function  of  the  parameter  i/ 
in  Fig.  IV. 4. 9.  It  appears  from  this  figure  that,  for  all  practical  purposes, 
the  Chapman-Jouguet  state  can  be  considered  as  well  established  by  the  time 
J  =  1.  This  means  that  the  Chapman-Jouguet  condition  is  practically  achieved 
when  y}\  =  .  This  radius,  according  to  our  definition  of  ,  is  deter¬ 

mined  by  the  value  of  the  initiation  energy,  £•  ,  rather  than  by  the  exother¬ 
mic  energy  of  the  detonating  medium,  .  In  view  of  the  fact  that  the  amount 
of  overdrivenness  of  the  detonation  wave  is  primarily  affected  by  the  initia¬ 
tion  energy,  this  result  is  mot  altogether  surprising. 

The  results  which  are  obtained  here  for  the  inviscid  solution  can  be  compared 
with  those  obtained  by  Kamel  et  al.  (1979).  The  same  problem  is  solved  by  two 
different  methods.  Kamel  used  the  phase  space  method  of  solution,  while  in 
this  worK  the  quasi-similar  technique  of  Oshima  is  used.  The  results  in  the 
two  cases  are  qualitatively  similar.  The  quasi-similar  technique  fails  to  ob¬ 
tain  the  values  of  the  gasdynamic  parameters  at  the  Chapman-Jouguet  state  which 
represents  a  singular  point  in  the  solution  domain.  Therefore,  the  solution  is 
stopped  at  =  0.033  which  is  approaching  the  value  of  0.033786.  However, 
the  quasi-similar  technique  is  easier  to  perform  and  it  is  faster  to  yield  re¬ 
sults  m  a  general  qualitative  manner. 

Figures  from  IV. 4, 10a  to  IV. 4. 17  show  the  comparison  between  the  results  of  the 


two  methods  for  the  pressure  profiles,  the  density  profiles,  the  temperature 
profiles,  the  velocity  profiles,  the  integral  curves  in  the  f  -  Z  phase 
plane,  the  energy  integral,  ^  ,  the  decay  coefficient,  A  *  and  the  non- 
dimensional  radius,  J  ,  respectively,  for  different  values  of  $  .  From 
these  figures,  the  two  methods  of  solution  appear  close  enough  to  warrant  the 
use  of  the  quasi-similar  approximation. 
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Figure  Captions 


Fig.  IV. 4.1  (a),  (b)  S  (c) 

Non-dimensional  pressure  profiles  for  planar,  cylindrical  and 
spherical  waves,  respectively,  at  different  values  of  y  with 
\  =  20.0231  and  S  -  1.3. 

Fig.  IV. 4. 2  (a),  (b)  &  (c) 

Non-dimensional  density  profiles  for  planar,  cylindrical  and 
spherical  waves,  respectively,  at  different  values  of  with 
\  =  20.0231  and  S'  =  1.3. 


Fig.  IV. 4. 3  (a],  (b)  &  (c) 

Non-dimensional  temperature  profiles  for  planar,  cylindrical 
and  spherical  waves,  respectively,  at  different  values  of 
with  =  20.0231  and  S'  =  1.3. 


Fig.  IV. 4. 4  (a).,  (b)  $  (c) 

Non-dimensional  velocity  profiles  for  planar,  cylindrical  and 
spherical  waves,  respectively,  at  different  values  of  ^  with 
\  =  20.0231  and  =  1.3. 

Fig.  IV. 4. 5  (a),  (b)  &  (c) 

Integral  curves  in  the  F-Z  phase  plane  for  planar,  cylindrical 
and  spherical  waves,  respectively,  at  different  values  of  ij 
with  =  2G.0231  and  =  1.3. 

Fig.  IV. 4. 6  Ness  integral  3*  as  a  function  of  with  =  20.0231  and 
S’  a  1.3  while  j  =  0,  1  and  2. 


Fig.  IV. 4. 7 


Energy  integral  as  a  function  of  ij  with  ^  =  20.0231  and 
S’  =  1.3  while  i  =  0,  1  and  2. 
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Fig.  IV. 4. 8  Decay  coefficient  as  a  function  of  ^  with  ^  =  20.0231  and 
X  =1.3  while  i  =  0,  1  and  2. 

Fig.  IV. 4. 9  Non-dimensional  radius  J  as  a  function  of  i j  with  %  =  20.0231 
and  $  =  1.3  while  J  =  0,  1  and  2. 

Fig.  IV. 4. 10  (a],  (b)  &  (c) 

Comparison  between  the  Quasi-Similar  and  Phase-Spac.e  methods  of 
solution  for  the  non-dimensional  pressure  profiles  with 

=  20.0231  and  X  =  1.3  while  j  -  0,  1  and  2,  respectively. 

Fig.  IV. 4. 11  (a),  (b)  &  (c) 

Comparison  between  the  Quasi-Similar  and  Phase-Space  methods  of 
solution  for  the  non-dimensional  density  profiles  with 

-  20.0231  and  X  =  1.3  while  d  =0,1  and  2,  respectively. 

Fig.  IV. 4. 12  fa].  (bJ  &  Cc] 

Comparison  between  the  Quasi-Similar  and  Phase-Space  methods  of 
solution  for  the  non-dimensional  temperature  profiles  with 
^  =  20.0231  and  X  =  1.3  while  d  =  0,  1  and  2,  respectively. 

Fig.  IV. 4. 13  (a),  (b)  &  (c) 

Comparison  between  the  Quasi-Similar  and  Phase-Space  methods  of 
solution  for  the  non-dimensional  velocity  profiles  with 
%  =  20.0231  and  X  =1.3  while  d  =  0,  1  and  2,  respectively. 

#, 

Fig.  IV. 4. 14  (a),  (b)  S  (c) 

Comparison  between  Quasi-Similar  and  Phase  Space  methods  of  so¬ 
lution  for  the  integral  curves  in  F-Z  phase  plane  with 
^  =  20.0231  and  =  1.3  while  i  =  0,  1  and  2,  respectively. 

Fig.  IV. 4. 15  Comparison  between  Quasi-Similar  and  Phase  Space  methods  of  so¬ 
lution  for  the  energy  integral  3^  with  %  =  20.0231  and 
X  =  1.3  while  j  =0,1  and  2. 
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Fig.  IV. 4. 16  Comparison  between  Quasi-Similar  and  Phase  Space  methods  of  so¬ 
lution  for  the  decay  coefficient  with  \  =  20.0231  and 
=1.3  while  4  =  0,  1  and  2. 

Fig.  IV. 4. 17  Comparison  between  Quasi-Similar  and  Phase  Space  methods  of  so¬ 
lution  for  the  non-dimensional  radius  'J  with  ^  =  20.0231  and 
S’  =  1.3  while  4  =  0,  1  and  2. 
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IV. 5.  ANALYSIS  OF  REACTIVE  BLAST  WAVES  PROPAGATING  THROUGH 
GASEOUS  MIXTURES  WITH  A  SPATIALLY  VARYING  HEAT  OF 
DETONATION* 

i)  Introduction 

In  the  actual  situation,  as  occurring  during  an  accidental  spill  of  fuel  gases 
into  open  air,  the  initial  conditions  are  no  longer  uniform  and  may  be  very  com¬ 
plex.  In  particular,  the  concentration  of  the  fuel  gases  may  vary  with  the  di¬ 
stance  from  the  fuel  source.  This  non-uniformity  will  affect  the  propagation  of 
any  resulting  detonation  wave  considerably  and  is  essential  in  studying  real  ex¬ 
plosions,  Geiger  (1979). 

In  this  section,  the  propagation  of  reactive  blast  waves  through  a  medium  of 
spatially  non-uniform  fuel  concentration  is  analyzed.  For  simplicity,  the  fuel 
gas  concentration  distribution  is  assumed  to  be  symmetrical  with  respect  to  the 
center  of  explosion.  With  this  assumption,  blast  wave  propagation  is  a  spatial¬ 
ly  one-dimensional  problem.  In  addition,  the  change  in  concentration  is  expres¬ 
sed  by  the  change  in  energy  released  at  the  wave  front.  The  chemical  heat  re¬ 
lease  ^  is  determined  by  the  chemical  composition  and  the  thermodynamic  states 
of  the  reactants  and  products.  In  this  analysis  %  is  assumed  to  depend  only  on 
the  initial  chemical  composition  which  varies  with  distance  from  the  explosion 
center.  While  the  transient  diffusion  problem  may  be  used  to  determine  the  con¬ 
centration  profile,  the  distribution  \  (  \r  ),  with  r  being  the  distance  from 
the  fuel  source  is  given  a  priori.  For  simplicity,  changes  in  concentration  due 
to  transient  diffusion  are  neglected  since  the  diffusional  velocities  are  much 
smaller  than  the  shock  wave  velocity. 

The  effect  of  combustion  heat  release  appears  only  in  the  boundary  conditions  at 
the  wave  front  and  varies  as  the  wave  propagates,  so  that  the  boundary  condi¬ 
tions  include  the  unknown  function  which  controls  the  propagation  of  the  wave. 
The  formulation  of  the  problem  is  applicable  for  the  three  geometrical  symme¬ 
tries:’  spherical,  cylindrical  and  planar,  while  the  results  are  obtained  in  the 
case  of  spherical  waves  only.  The  detonating  medium  is  assumed  to  behave  as  an 
invibcid  perfect  gas  with  a  mean  specific  heat  ratio  S’ 


This  application  is  based  on  Ohyagi  et  al.  (1961). 
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ii)  Problem  Formulation 

At  a  time  t  =  0  and  a  point  V"  =  0  a  finite  amount  of  energy  is  liberated 
by  some  energy  source  to  form  a  non-steady,  symmetrically  expanding  shocK  wave, 
i.e.  a  blast  wave  as  is  shown  in  Fig.  IV. 5.1.  If  the  medium  is  an  inert  gas, 
the  blast  wave  decays  to  sound  wave.  In  exothermic  reactive  media,  if  the 
energy  exceeds  a  critical  level,  the  waves  decay  to  the  Chapman-Jouguet 
wave  whose  propagation  velocity  is  determined  by  the  chemical  energy  of  the  me¬ 
dium. 

In  this  problem,  the  concentration  of  the  fuel  gas  changes  with  the  distance 
from  the  center  of  explosion.  As  noted  previously,  the  distribution  utilized 
here  is  that  of  chemical  energy  instead  of  that  of  the  concentration.  For  sim- 

I 

plicity,  the  distribution  is  assumed  to  be  symmetric  with  respect  to  the  center 
so  that  a  symmetrical  nature  of  the  propagation  is  assumed.  The  chemical  ener¬ 
gy  release  at  the  wave  front  is  given  as  a  function  of  i£  , 

%  -  <IV*5*1J 

Any  explosive  gas  mixture  exhibits  detonability  limits  beyond  which  a  detona¬ 
tion  can  not  be  established.  Here,  the  chemical  reaction  is  always  assumed  to 
be  complete  at  the  wave  front,  and  the  complex  physico-chemical  phenomena  oc¬ 
curring  in  the  reaction  zone  is  neglected.  If  the  detonatility  limits  should 
be  taken  into  account,  the  value  of  the  limiting  concentration  or  the  limiting 
chemical  energy  should  be  given  a  priori. 

In  addition  to  these  fundamental  assumptions,  one  may  assume  that  the  change 
in  molecular  weights  of  the  fuels  and  the  oxidizers  do  not  affect  the  density 
and  the  pressure  of  the  medium  ahead  of  the  front. 

The  governing  equations,  Eqs.  (11.20),  (11.22)  and  (11.24),  for  a  sourcelsss, 
:;>viscid  and  thermally  and  calorically  perfect  with  constant  specific  heats 
medium,  are  reduced  to 
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2£  +rA»  +  X  =o  L  (iv.5.2) 

M  $|  +(/-_*>2i-A9+^e£*i|)=o  _ 

where  X  is  the  decay  coefficient  defined  in  Eq.  (11.14),  which  may  be  rewritten 
as 


^  {4  -  _  2  ^  1*  H 

diVufi  cJ&ift 


(IV. 5. 3) 


The  boundary  conditions  of  the  problem  are  that  of  the  reactive  medium  with  in¬ 
stantaneous  hBat  release  at  the  wave  front  of  strong  detonation,  which  are  given 
by  Eqs.  (11.93),  (11.54)  and  (11.85).  These  equations  may  be  rewritten  as  fol¬ 
lows  : 


where 


f  ,  fa,  •d)  . 
n  *>1 


kn  =  h (i,V)  «  —£±i  - 

9<i,U)  3  1 


(IV. 5. 4) 


4  V« 

K  =  £(i-  tf)  -  2 


(IV. 5. 5) 


with  ^  being  the  non-dimensional  chemical  heat  release  at  the  wave  front,  de- 
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fined  by  Eq.  (11.53)  which  is,  in  this  case, 

?  =  (IV. 5. 6) 

< 

The  heat  release  at  the  wave  front,  %  ,  depends  on  the  wave  position,  , 
so  that  it  is  a  function  of  . 

A  condition  to  be  satisfied  at  the  center  of  symmetry  is  that  the  particle  ve¬ 
locity  must  equal  to  zero  for  all  times,  which  is  given  by  Eq.  (11.112)  as 

f(Oj  i 1)  ~  0  (IV. 5. 7) 


The  problem  now  is  to  solve  Eqs.  (IV. 5. 2)  with  the  boundary  conditions,  Eqs. 
(IV. 5. 4),  as  well  as  the  compatibility  equation,  Eq.  (IV. 5. 7),  to  find  a  func¬ 
tional  form  of  for  a  given  function  of  .  For  this  purpose 

the  relation  between  rn  and  must  be  given. 

Integration  of  Eq.  (IV. 5. 3)  yields 


y 

}  s  ir  s  e*P  (  (  )  (IV. 5.8) 

where  d0  and  ^  are  arbitrary  constants  which  express  the  initial  condition 
of  the  wave  Mach  number  and  position,  respectively.  When  the  value  of  *d0  is 
very  small,  the  self-similarity  of  the  flow  holds  at  this  initial  stage.  As 
defined  previously  in  Eq.  (11.72),  the  initial  position  l£  is  proportional 
tQ  Bt  ,  By  using  Eq.  (IV. 5. 3),  may  be  expressed  in  terms  of  ^  . 

Therefore,  the  boundary  conditions,  through  cj.  )  ,  as  well  as  the  basic 
conditions  include  the  unknown  function 

here,  the  function  %  IK)  is  assumed  to  be  differentiable.  Therefore, 
the  application  of  the  quasi-similar  concept  yields  the  following  ordinary 
differential  equations  for  the  governing  equations,  as  given  by  Eqs.  (III. 23)- 


(III. 25)  where,  for  a  sourceless  flow  field,  all  ^  ’s  vanish. 


(A  x)  ik.  +  k(d£  +  i£)  +  %f\  h  —  o 

d  X'  oX 


vhile 


(f-x)M.  +±-d I  +%&f  = 

Jx  h  6X 


(IV. 5. 9) 


where 


| +  **$+*£>+  ^ c  j 


(IV. 5. 10) 


8  s  .5?L  —  0.5 


(IV. 5. 11) 


y<^  f-  ^  _  1 


(IV. 5. 12) 


a  4#-  =*^4>  (44  -*> 


(IV. 5. 13) 


Plll«l 
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(IV. 5. 14) 


As  an  example,  a  simple  exponential  form  of  the  function  C»ft)  is  adapted. 
In  the  non-dimensional  form,  it  is 


exp  (-p  J&.) 


exp  (-p  J  ) 


(IV. 5. 15) 


where  ^  and  /3  are  non-dimensional  parameters. 

If  /3  is  positive,  the  fuel  concentration  becomes  smaller  with  the  distance  from 
the  center  and  the  chemical  energy  released  at  the  wave  front  decreases  exponen¬ 
tially  as  the  wave  propagates. 

In  the  case  of  Chapman-Jouguet  detonation,  K  equals  to  zero,  see  Eq,  (IV. 5. 5), 
thus 


Cf.  -  (*  - 


To  obtain  a  functional  form  for  the  decay  coefficient  in  the  case  of  C-J  wave, 
one  may  differentiate  the  above  equation  with  respect  to  'J  and  utilize  the  de¬ 
finition  of  ^  ,  given  by  Eq.  (IV. 5. 3),  to  get 


A  =  ,<Lkl)  ,  ?(*-*>% 

CJ  Jf'ij'cj  Jlfi-'i)  +(»*-«?] 


(IV. 5. 16) 


iii)  Solution 


The  set  of  simultaneous  ordinary  differential  equations,  Eqs.  (IV. b. 9)  is  solved 


numerically.  Since  AW  is  not  Known  initially,  the  value  of  A  for  each  y 
is  first  assumed  and  integrations  are  started  from  the  wave  front,  9C  =  1, 
where  the  boundary  conditions  are  given  by  Eqs.  (IV. 5. 4)  to  the  center  of  ex¬ 
plosion,  %  -  0,  where  the  particle  velocity  is  equal  to  zero.  If  this  condi¬ 
tion  is  not  satisfied  for  the  assumed  value  of  A  »  this  value  is  corrected 
and  the  procedure  is  repeated  until  Eq.  (IV, 5. 7)  is  satisfied  to  the  required 
accuracy. 

A  difficulty,  however,  arises  because  the  integrated  term  of  A(tf)  appears  in 
the  boundary  conditions.  The  steps  of  the  solution  in  this  case  are  as  fol¬ 
lows  : 

a)  To  obtain  Aod)  ,  y  is  assumed  firstly  to  be  equal  to  . 

bJ  From  Eq.  (IV. 5. 8],  J  =  1  and  is  calculated  from  Eq.  (IV. 5. 6)  so  that 
the  decay  coefficient  Xtfo)  is  obtained  by  iteration. 

c)  Next,  y  is  increased  by  small  amount,  A  y  ,  and  ^  is  calculated  by  nu¬ 
merical  integration  of  Eq.  (IV. 5. 8)  with  the  trapezoidal  rule. 

d)  Applying  the  same  procedure  as  described  |i,  Piously  leads  to  the  solu¬ 
tion  A  0*o  +  Ay)  . 

d)  The  further  solution  will  be  found  by  lapel ition  of  this  procedure. 

It  should  be  noted  here  that  the  error  in  this  procedure  can  be  reduced  if  the 
increment  is  small. 

iv)  Results  and  Conclusions 

The  results  obtained  here  are  for  the  case  of  spherical  wave,  J  -  2.  The  va¬ 
lues  of  the  other  parameters  used  are: 

=  20,  p  =  0.001,  0.05,  0.1,  0.15,  0.2  and  0.5,  y„  «  0.0001  and  2T  =  1.4. 

So,  in  these  calculations,  the  medium  has  suf t iciont,  chemical  energy  to  support 
a  Jetonation  near  the  explosion  center,  with  u,..,  energy  released  decreasing  ex¬ 
ponentially  with  distance. 

It  should  be  noted  here  that  if  /9  is  equal  to  iLiu,  the  wave  will  decay  to  the 
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Chapman- Jouguet  wave  with  =  0.026,  calculated  by  using  Eq.  [11.94),  and 

the  problem  is  reduced  to  that  discussed  briefly  in  section  (IV. 4). 

i 

The  step  size  for  the  Runge-Kutta  integration  for  %  -  direction,  A  X  ,  is 
chosen  to  be  0.01.  For  the  y  -  direction,  the  step  size  A^  is  0.0001. 

Figure  IV. 5. 2  shows  the  decay  coefficient  X('d)  for  each  J3  .  When  y  is 
equal  to  zero,  the  self-similar  solution  is  fulfilled  so  that  X  is  equal  to 
j+  i  =  3.  As  y  increases,  X  decreases  generally  which  means  that  the 
wave  approaches  the  steady  wave.  For  small  values  of  p  ,  as  0.001,  the  heat 
liberated  at  the  wave  front  is  almost  constant  and  the  wave  behaves  as  a  re¬ 
active  blast  wave  in  the  uniform  combustible  mixture  for  the  regime  where  the 
solution  can  be  found.  It  seems  to  decay  to  the  C-J  wave  (  =  0.026)  but 

the  solution  can  not  be  found  near  C-J  state  because  of  the  singularity  at 
the  wave  front  in  the  C-J  wave.  In  the  C-J  wave,  the  gradients  of  p  ,  k 
and  9  in  the  -  direction  become  infinity  as  indicated  previously  in  auc¬ 
tion  (IV. 4).  As  p  increases,  l#l  decreases]  and  a  transition  from  the  re¬ 
active  blast  wave  to  the  non-reactive  one  occurs  with  increasing  9  .  For 
/3  =  0.05,  X  has  a  minimal  value  of  1.09  at  the  value  of  ^  a  0.1)5.  It 
appears  that  the  wove,  decaying  as  a  reactive  blast  wave,  loses  its  energy  as 
it  propagates  and  in  turn  decays  more  rapidly  as  a  non-reactive  blast  wave. 

But  for  this  case,  the  solution  can  not  be  found  for  3  >■  0.07  because  of  the 
limit  of  the  numerical  precision.  For  p  =  0.1,  0.15,  0.2  and  0.5,  the  com¬ 
plete  solutions  show  the  transition  from  reactive  blasts  to  non-reactive  ones. 
For  p>  =  0.5,  the  solution  is  almost  identical  to  that  for  ^  -  0  (or  p  -»«»), 
i.e.,  the  non-reactive  blast  wave  except  for  very  small  values  of  y  .  The 
values  of  the  decay  coefficient  for  /S  >  1,  when  9  >.  0.35,  or  M  <  1.7,  are 

nearly  the  same  because  the  chemical  energy  2j.  ,  there,  is  almost,  equal  to 
zero.  The  dashed  curve  shown  in  this  figure  represents  the  Chapinan-Jouguet  de¬ 
cay  coefficient,  •  for  J3  =  0.1,  as  a  function  of  9  which  is  given  by 

Lq.  (IV.5.16J.  The  C-J  decay  coefficient  for  f3  =  0.1  varies  with  9  because 
the  chemical  energy  ^  changes  with  J  ,  and  it  can  be  called  as  the  local  C-J 
decay  coefficient. 

Figure  IV. 5. 3  shows  that  the  shock  l"ach  number,  M  ,  decreases  with  an  increase 
in  ih-_  non-dimensional  radius,  T  .  The  curves  for  all  1  s  are  calculated 
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with  the  same  initial  condition,  i.e.,  for  J  =  1,  l“1  =  100.  For  p  =  0.001, 
the  Mach  number  ti  seems  to  decay  to  the  Chapman-Jouguet  value,  =  6.202, 

corresponding  to  ^  qj  =  0.026.  In  reality,  it  decays  to  the  sound  wave  as 
the  wave  propagates  infinitely  apart  from  the  center.  As  p  increases,  f'l 
seems  to  decrease  to  unity  in  the  far  field.  From  Fig.  IV. 5. 3,  it  can  be  seen 
that  for  "J  >.  25,  fi  is  less  than  1.7  for  >  0.1.  Therefore,  in  the  re¬ 

gime  J  25,  although  q.  is  almost  zero  for  those  /3  ’s,  ^  for  larger  P 
is  greater  than  that  fur  smaller  {3  's.  It  can  be  said  that  this  is  the  effect 
of  history  which  the  wave  has  experienced.  The  dashed  lines  in  Fig.  IV. 5. 3  in¬ 
dicate  the  Chapman-Jouguet  Mach  number  obtained  by  K  =0.  The  C-J  Mach  number 
for  a  certain  p  varies  with  J  because  \  is  changing  and  it  can  be  considered 
as  the  local  C-J  Mach  number.  It  is  evident  that  the  Mach  numbers  for  all  p  's 
decay  to  their  own  local  C-J  values. 

Figures  IV. 5. 4  through  IV. 5. 6  show  the  variations  of  the  peak  values,  just  be¬ 
hind  the  wave  front,  of  the  non-dimensional  density,  velocity  and  pressure  with 
respect  to  J  .  They  change  from  the  values  of  the  strong  blast  wave  limits 
(i.e.  I -  6.0,  -  0.8333  and  =  0,8333  for  K  =  1.4)  to  the  values  of 

the  C-J  wave  (i.e.  !%  =  1.683,  s  0.4508  and  =  0.4244  for  =  1.4  and 

=  20)  for  the  reactiVe  blast  wave  in  the  uniform  medium,  while  they  ap¬ 
proach  the  values  of  sound  waves  (i.e.  Hrt  =  1.0,  fn  -  0  and  =  ^0  for  the 

wave  in  the  non-uniform  medium.  From  Fig.  IV. 5. 4,  it  can  be  said  that  in  the 
early  stage  the  peak  density  at  the  wave  front  decays  more  rapidly  for  small  p 
than  for  large  (3  with  respect  to  the  distance  from  the  center.  In  the  early 
stage,  M  is  determined  by  the  initiation  energy  rather  than  the  chemical  ener¬ 
gy.  For  the  same  M  ,  the  peak  values  for  small  chemical  energy  is  larger 
than  for  large  q.  as  is  calculated  by  the  Hugoniot  relation.  In  the  later 
stage  of  propagation,  approaches  to  the  local  C-J  value  and  it  is  determined 
by  the  chemical  energy.  Therefore,  the  peak  density  for  large  p  has  a  tenden¬ 
cy  to  decay  more  rapidly  than  that  for  small  p  . 

figures  IV'.5.7a-c  show  the  typical  density  profiles  normalized  with  respect  to 
the  peak  values  at  the  wave  front  for  p  oo  (or  ^  =  0),  fi  =  0.2  and  0.001. 
Tnu  non-dimensional  density  near  the  center  increases  as  ^  1  for  p  so  and 

P  =  j.2  (here,  the  profiles  of  p  -  0.2  for  $  >  0.25  are  omitted  because 


they  are  almost  the  same  as  those  for  p  <»  ).  For  P  =  0.001,  the  non- 
dimensional  density  increase  very  little  as  ^  (when  (d  =  0.022, 

=  0.646).  In  these  calculations  of  the  point  blast  wave,  the  density 
near  the  center  remains  zero,  corresponding  to  the  fact  that  the  velocity 
gradient  a£.  does  not  vanish.  If  it  is  possible  to  extend  the  solution  for 

T)  % 

larger  ^  ,  the  density  may  increase  more  near  the  center. 

Figures  IV.5.8a-c  show  the  particle  velocity  profiles.  The  core  in  which 
the  velocity  vanishes  is  not  accurate  in  these  figures.  This  is  a  consequence 
of  the  quasi-similarity  approximation  since  it  is  expected  that  for  -  0 
and  >  0.6  then  f  w  has  a  negative  region  near  the  center. 

Figures  IV.5.9a-c  show  the  pressure  profiles  for  p  -*r  oo  ,  /3  =  0.2  and 

^  =  0.001,  respectively,  where  the  piofiles  exhibit  the  same  general  trend 

( 

as  in  previous  cases. 

In  this  model  using  the  quasi-similar  approximation,  there  have  to  be  some  in¬ 
accuracies  in  the  gasdynamic  profiles.  However,  the  decay  coefficient  obtai¬ 
ned  seems  to  be  reasonable.  In  addition,  this  method  permits  the  use  of  any 
model  for  the  heat  release  function  and  it  makes  it  very  convenient 

to  find  the  decay  coefficient  for  the  wave  front. 

In  this  example,  is  taken  to  decrease  exponentially  with  distance  from  the 
center,  the  nature  of  wave  propagation  is  determined  by  the  ratio  of  the  cha¬ 
racteristic  radius  of  the  initiation  energy,  lr£  ,  to  that  of  the  decreasing 
initial  chemical  energy,  J3  ,  and  it  is  revealed  that  for  small  ti  «  such  as 
P.00:  and  2.C5,  the  wave  front  behaves  as  a  detonation  wave  and  for  large  va¬ 
lues  of  p  ,  such  ,ns  0.5,  it  decays  as  a  non-reactive  blast  wave. 


Fig.  IV. 5.1  Configuration  of  a  spherical  blast  wave. 

Fig.  IV. 5. 2  The  decay  coefficient  H  as  a  function  of  if  for  different  values 
of  p  with  If  =  1.4  and  ^  =  20  while  j  =  2. 

Fig..  IV. 5. 3  The  shock  Mach  number  M  as  a  function  of  the  non-dimensional 
shock  radius  J  for  different  values  of  p  with  if  =1.4  and 
*1  =  20  while  j  *  2. 

Fig.  IV. 5. 4  The  density  ratio  across  the  wave  front  as  a  function  of  the  non- 
dimensional  shock  radius  J  for  different  values  of  p  with 
=1.4  and  <L  ■  20  while  i  ■  2. 

Fig.  IV. 5. 5  The  velocity'  ratio  across  the  wave  front  as  a  function  of  the  non 
dimensional  shock  radius  J  for  different  values  of  p  with 
S'  *  1.4  and  4.  =  2!J  while  j  =  2. 

Fig.  IV. 5. 6  The  non-dimensional  pressure  at  the  wave  front  Sn  as  a  function 

of  the  non-dimensional  shock  radius  J  for  different  values  of  p 
with  ^  =1.4  and  ^  =  20  while  j  =2. 

Fig.  IV. 5. 7  (a),  (b)  S  ,c) 

The  density  profiles  at  different  values  of  ^  for  p  =  oo  ,  3.2 

and  O.QO'I,  respectively,  with  If  =  1.4  ana  4  =20  while  J  =  2. 

ro 

t  1.  .  iv.it.d  fa),  (bj  &  fc) 

The  particle  velocity  piofiles  at  different  values  of  y  for 
P  =  <=o  ,  0.2  and  0.001,  respectively,  with  If#  1.4  and  ^  =  20 

while  J  -  2. 

r  ....  I.  .  ..1  ClO  &  Ic) 

fhe  pressure  profiles  at  different  values  of  $  for  P  -  , 


0,2  and  0.001,  respectively,  with  if  =  1.4  and  ^  =  20  while 

’© 

6  =  2. 
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IV. 6  BLAST  WAVES  IN  A  DETONATING  MEDIUM  WITH  TRANSPORT 
PROPERTIES  TAKEN  INTO  ACCOUNT* 

i)  Introduction 

It  is  evident  from  previous  applications,  for  both  reactive  and  non-reactive 
blast  waves,  that  the  solution  fails  near  the  center  of  symmetry  since  the 
temperature  goes  to  infinity,  while  the  density  of  the  gas  tends  to  zero.  In 
order  to  improve  the  solution  neai  the  center,  one  has  to  taka  the  effects  of 
viscosity, heat  conduction  and  radiation  into  account. 

In  this  section,  the  non-self-similar  problem  of  a  point  explosion  in  a  deto¬ 
nating  gas  with  viscosity,  thermal  conduction  and  radiation  effects  is  analy¬ 
zed.  Radiation  has  been  included  by  assuming  a  diffusion  radiation  mode  for 
an  optically  thick  grey  gas  which,  as  a  consequence,  leads  to  radiation  terms 
that  are  mathematically  similar  to  those  of  thermal  conduction.  The  trans¬ 
port  coefficients  are  assumed  to  be  proportional  tc  appropriate  powers  of  tem¬ 
perature.  A  boundary  layer-like  region  is  assumed  to  exist  near  the  blast 
origin  where  all  the  transport  effects  are  concentrated,  thus  permitting  the 
rest  of  the  flow  field  to  be  treatea  in  the  usual  inviscid  manner.  After  the 
inviscid  solution  has  been  obtained,  series  expansions  of  the  gasdynamic  para¬ 
meters  are  then  applied  to  the  inner  region  equations  which,  in  turn,  reduce 
to  algebraic  relations.  The  structure  of  the  non-self -similar  flow  fields  is 
then  fully  determined  by  matching  the  approximate  solution  for  the  inner  re¬ 
gion  with  that  for  the  inviscid  outer  region.  The  three  geometrical  symme¬ 
tries  of  blast  waves  are  considered  and  the  detonating  medium  will  be  treated 
as  a  perfect  gas  with  a  mean  specific  heat  ratio,  . 

^ )  Problem  Formulation 

The  three  conservation  equations  of  mass,  momentum  and  energy,  taking  into  ac¬ 
count  possible  sources  of  mass,  momentum  or  energy  that  may  affect  the  flow 
field,  are  given  by  Eqs.  (II. 1)  -  (II. 3). 

For  me  case  considered  here,  the  rate  of  mass  supplied  per  unit  mass  of  the 
mixture  is  zero.  Thus 


*  Inis  application  is  basec  on  Abdel-Raouf  (1982). 
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(IV. 6. 13 


when  taking  the  effect  of  thermal  radiation,  by  assuming  a  diffusion  radiation 
mode  for  an  optically  thick  grey  gas,  thermal  conduction  and  viscous  forces 
into  account,  the  momentum  and  energy  source  terms,  as  given  by  Kamel  et  al. 
(1977a)  take  the  form: 


r1 


+' i(^  _  W.)l  +  JU 

r  '  r  ' J  3  f  %r  ' \r  -4  r 


(IV. 6. 2) 


and 


yi 


+  Jl  3X1  (U  _  X_  Ux)  +  Ji.  (r*  3X  ) 

3  ^  r  '  ir  /  r’j  r-  '  ~s>r  ' 


^  »t“3 


+  XU  .  Xi  +  JU  (JAJLJL  )  X_  /r d  XL  ) 

r  ~sr  r  J  v  3  CL  '  *° r  V  r 


/  16  s'  T*  \  ~*T 


/  /  6  XX  ) 

-Tv  ^  ^  >  * 

3  a 


~t>  r- 
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where  T  is  the  absolute  temperature  and  <r'  is  the  Stefan-Boltzmann  constant, 

✓ 

while  JJ-  ,  K  and  O-  are  the  coefficients  of  viscosity,  thermal  conduction  and 
absorption,  respectively.  These  coefficients,  being  thermodynamic  properties, 
are  functions  of  temperature  and  pressure.  However,  the  pressure  is  very  near¬ 
ly  constant  in  the  inner  region  near  the  center  of  symmetry  where  the  transport 
properties  predominate.  Therefore,  these  coefficients  may  be  considered  as 
functions  of  temperature  only  and  could  be  assumed  to  vary  with  the  local  tem¬ 
perature,  as  given  by  Kamel  et  al. (1973a),  in  the  following  manner: 


/*  */*oO. f)  j  K  -  Ko<£> 


a  (-&  ) 
u 


(IV. 6. 41 


where  cn  and  p  are  constants,  and  subscript  o  refers  to  some  Known  state. 

One  may  then  define  the  following  non-dimensional  parameters  of  the  transport 
properties 


,  Brm  /U°C± 


Ho  je 


4  c r  T1 

o 


(IV. 6. 5) 


wiiei  e  R.k  is  the  Reynolds  number;  ,  the  Prandtl  number;  No  ,  the  radia¬ 
tion-conduction  parameter;  Cp  the  gas  specific  heat  at  constant  pressure  and 
to  is  the  characteristic  length  of  the  explosion. 

The  basic  equations  of  mass,  momentum  and  energy,  in  non-dimensional  form,  Eqs. 
(11.20 ) ,  (II. 2?)  and  (11.83),  with  the  definitions  of  Eqs.  (11.21),  (11.23) 
and  (11.84),  under  the  assumption  of  uniform  ambient  atmosphere,  are  then  re¬ 
duced  to 


M  .  jit  +  4  f  _ 


_  A  +  M  2l  (f-X)  ^  2d 

1  f  -*1  F  kf 
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where 


c±  ■  4  * 
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R*  J 


W  * 


4-P 
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R  #  iV  A/o  J 
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r«  *v  j 
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The  boundary  conditions  of  detonating  blast  wave  problems  are  given  by  the  Hugo- 
niot  relations,  Eqs.  (11.93),  (11.54)  and  (11.85).  These  equations  may  be  re¬ 
written  as 


1 


and 


(IV. 6. 11) 


j 


Trie  reduction  of  the  governing  equations,  Eqs.  (IV. 6. 6),  to  self-similarity 

leads  to  the  determination  of  the  values  of  the  constants  cc  and  pi  ,  as  given 

», 

by  Kanicl  ct  al.  (1877a),  which  aio.  equal  to;. 
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c< 


(IV. 6. 12) 


The  coefficients  of  viscosity,  thermal  conduction  and  absorption  which  vary 
with  temperature  as  given  by  Eq.  (IV. 6.4),  must  be  independent  of  time,  since 
they  are  properties  of  the  flowing  substance.  Therefore,  the  constants  oi.  and 
P  must  have  the  same  values,  given  by  Eq.  ( ik/ .6.12),  that  are  independent  of 
y  which  are  also  valid  for  the  non-self-similar  case. 

The  mass  and  energy  integrals,  in  the  case  of  inviscid  medium,  have  the  forms 
given  by  Eqs.  (11.65)  and  (11.96)  respectively,  which  are 


(IV. 6. 13) 


and 


+ 


1 

(Ui)  (8-1) 


(IV. 6. 14) 


while  the  decay  coefficient  X  is  given  by  Eq.  (11.97)  as 


d 

'X  „  <u  j)  $3  -  jit- 1 )  -  H-d 


or.  - 


a  d. 

*  diJ 


which  for  self-similar  cases,  reduces  to 


(IV. 6. 15 
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The  quasi-similar  technique,  when  applied,  reduces  the  basic  equations,  Eqs. 
(IV. 6. 6),  to  ordinary  differential  ones,  given  by  Eqs.  (III. 23)  -  (11.25). 
Noting  that  =  0,  these  equations  in  this  case  reduce  to 


(f-  x)  Ml  +  k  (J£  +  iJL)  +  'XhQ  =  O 
Ax  Ax  * 

US  =f*r 

and 


where 


and 


A  e 


<P 


o.  5 


Q.  «=  %  •+■  V  i 

with 

CP  =  y  ik  )  fi  . 

d^/  L  -(V+l)£,  J 

,  L  :  •.  n_  :ti -.utivec  ---  »  jik  ana 

6*  Ax 


i!  in  Eqs.  UV.C.l/l.  one  cb- 

A% 


tains,  as  given  by  Eqs.  (III. 26)  -  (III. 28)  when  -  0, 


A  =  r_L  (j£  +  j£  +  XA) 

A*  (f.X)  K  A*  *  x 


(IV. 6. 18) 


t 


and 

=  -  h  [cf-  *)J£  +  tfe  .  f<?  ] 

AX  ^  dX 


iii )  Solution 

The  solution  procedure  is  divided  into  two  stepB.  First,  the  inviscid  flow 
field  is  obtained  by  solving  the  conservation  equations,  Eqs.  (IV, 6. 18),  ne¬ 
glecting  ♦vansport  effects  (i.e.  <Pr  -  ■  0)  and  £he  problem  is  reduced 

to  the  one  discussed  previously  in  section  (IV. 4).  This  solution  is  obtained 
in  order  to  determine  X  .  Then  the  effect  of  transport  phenomena  is  accoun¬ 
ted  for  by  assuming  the  flow  field  tc  consist  of  two  regions,  namely  an  outer 
region  and  ar  inner  region,  with  transport  effects  existing  in  the  inner  re¬ 
gion  only  (Kamel  et  al.,  1977a).  Series  expansions  for  the  gasdynamic  parame- 
t  :  ^  aie  usua  to  obtain  the  inner  region  solution  which  be  matched  with 

i'  r  >r  the  inviscid  outer  region. 

/m,  .jiJui  ut  magnitude  analysis  to  the  conservation  equations  was  performed  by 
f  u:\el  wt  vi!*tll}7?a)  to  jimplify  them  in  the  inner  region..  Two  main  corclu- 

,  ■■  .  1  t. 

. )  1 nu  vijCd-ity  terms  are  of  negligible  contribution  and  may  thu-.  oe 
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dropped  from  the  momentum  and  energy  equation, 
b)  The  momentum  equation  rt duces  to: 


(IV. 6. 19) 


which  indicates  that  the  pressure  is  constant  throughout  the  inner  region.  The 
continuity  equation  remains  unchanged,  as  given  by  the  first  of  Eqs.  (IV. 6. 6), 
while  the  energy  equation  becomes: 


-  1 .  i-  +  (4  +  .*)  J_  \*  1 
*  h  tx  '  J 


(IV. 6. 20) 


where 


*  =  £i  +  CJ 

Fur  the  inner  region  solution,  instead  of  following  the  complicated  procedure 
of  finding  the  value  of  which  contains  second  order  differentiations  by 
iteration,  one  can  obtain  an  approximate  solution,  leading  to  at  least  a  good 
qualitative  description  of  the  flow  field,  by  assuming  that  the  gasdynamic 
parameters  may  be  expanded  in  a  power  series  in  the  form 


f  *  fo  +  +  o<x4 ) 


k  =  k0  +  /t,  X  +  kt  X1  +  Aa  X3+  0(X*)  |p 


(IV. tf. 21. 
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and 

9  =  9, 


where  the  coefficients  of  the  powers  of  X  are  functions  of  only.  A  com¬ 
patibility  equation  which  must  be  satisfied  by  the  correct  solution,  Eq. 
(11.112),  is  that  the  particle  velocity  at  the  center  of  symmetry  must  be 
equal  to  zero,  that  is  /(o.’ij)  =  0.  Immediately  it  is  shown  that  £  =  0. 

At  the  center  of  symmetry,  however,  one  has: 


-aT 


at 


X  a  O 


Equation  (IV. 6. 19)  with  the  above  equation  may  yield: 


2Jl 


at 


XtzO 


If  this  equation  is  to  be  satis- led  then  Wi  must  also  vanish.  Substituting 
Eqs.  (IV. 6. 21)  into  the  fiiat  of  Eqs.  (IV. 6. 6)  and  Eq.  (IV, 6. 20),  and  equa¬ 
ting  the  coefficients  of  the  same  powers  of  X  to  zero,  one  obtains  from  the 
continuity  equation: 

^  ^  T2  *  (IV. 6. 2: 

n.0 

v.r.-vVt  the  prime  indicates'  differentiation  with  respect  to  ,  also  F  s  o 
A  -•  o  unu 


X  d  (k'2  -  A  |  jig  ,  JhJ-  )  -  2  kt  +  (U  3)  A „P3=0 

rt  o  1 
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From  the  energy  equation,  one  gets 


,  ,  Q« 

+  |s.  -  -7=7,  (lv-s-24) 

90  ''-o  h* 

The  gasdynamic  parameters  /  and  A  .  given  by  E,a.  HV.6.21J.  are  now  reduced 
to: 


f  -  ft  X  *  £  x‘ 


and 


k 
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In  addition  to  algebraic  Eqs.  (IV. 6. 22)  -  (IV. 6. 24)  one  has  the  three  equa 
tions  that  match  the  inner  solution  to  the  outer  solution,  namely 


a  +  hx 

SL-  % 


(IV. 6. 26) 
P 


where  is  the  %  boundary  of  the  inner  region. 

Equations  (IV. 6. 22),  (IV. 6. 23),  (IV. 6. 24)  and  (IV. 6. 26)  form  a  system  of  six 
equations  that  may  be  solved  simultaneously  for  f  >  h0  •  k^  ‘ 

and  X  •  f°r  a  particular  shock  strength  !j  . 

fit 


The  derivatives  of  the  coefficients  k0  ,  hx  and  S0  with  respect  to  3  can 


be  put  in  a  finite  difference  form  as: 
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To  develop  *-he  computational  procedure,  the  problem  is  solved  first  in  the 
self-similar  case  by  putting  y  equal  to  zero  in  Eqs.  (IV. 6.22)  -  (IV. 6. 24). 
These  equations  reduce  to 


h  *  ° 


(IV. 6. 2Q) 


After  some  algebraic  manipulation,  one  can  easily  obtain  the  following  equa¬ 
tions  for  the  determination  of  the  coefficients  %  ■  F.  ■  k.  and  4*  in 
terms  of  the  values  of  the  gasdynamic  parameters  at  the  matching  paint  3^  , 
/  and  k  ■: 


£  *  L  / 

I 

W0~  z  /[(/+3)  /3  X^  +  z] 


and 

=  4  K  <J+3)  A 
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Since  all  the  coefficients  depend  ori  the  matching  point  parameters,  an  itera¬ 
tive  procedure  is  required  to  find  the  matching  point  which  satisfies  the 
last  of  Eqs.  (IV. 6. 20).  The  latter  equation  may  not  be  satisfied  at  a  speci¬ 
fied  step  of  the  numerical  integration  procedure,  JC(-  .  Thus,  assuming  a  li¬ 
near  relation  between  the  residual  value  of  that  equation,  R  ,  and  X  ,  the 
correct  value  of  the  matching  point  can  be  obtained  by  interpolation,  as 
shown  in  Fig.  IV.6.1.  In  addition,  the  following  relations  were  used: 


*  -  a,  -  C‘/(‘<  C) 


(IV.6.3U: 


and 


*n  =  R,7di(-  ft,.,  ) 


(IV. 6. 31) 


where  R(  is  the  residual  value  of  the  last  of  Eqs.  (IV. 6. 28)  at  %■  . 

The  corresponding  values  of  the  gasaynamic  parameters  at  the  matching  point 
are  also  determined  by  interpolation,  considering  straight  line  relations  bet¬ 
ween  their  values  and  %  .  For  an  example.  Fig.  IV. 6. 2  shows  how  the  particle 
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velocity  at  the  matching  point,  f  ,  is  determined.  Thus 


f  =  C  +  (fi.i  -P.-XXn-*^ 


k  _  +  ~ki)(X„-Xt) 

ftt  $  '  ..^r.v-r-n  , . 

(AX) 


(IV. 6. 32} 


<9  _  q  (9,U  -  9i)(Xm  -M 

(AX) 


Once  the  matching  point  parameters  are  determined,  the  coefficients  %  •  k.0  • 
ki  ana  are  also  determined  to  define  completely  the  inner  region  solu¬ 
tion  for  me  self- similar  problem. 

One  can  then  proceed  to  obtain  the  solution  for  the  non-self-simllar  problem. 

At  a  specified  value  of  d  =  ,  from  Eqs.  (IV. 6. 221,  (IV. 6. 23)  and  (IV. 8. 27), 

after  some  algebraic  manipulations,  one  may  obtain  the  coefficients  %  ,  |x9  , 

ki  ,  fx  and  f3  in  terms  of  the  values  of  the  gasdynamic  parameters  at  the 
matching  point  as  follows: 
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where 


=  Ifn-ft,  *»«)  /K 
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Equation  (IV. 6. 24)  must  be  satisfied  for  the  correct  values  of  the  gasdynamic 
parameters,  and  thus  one  follows  the  same  procedure  which  is  used  for  the  self¬ 
similar  problem. 

The  equation  which  corresponds  to  Eq.  (IV. 6. 30)  will  take  the  form: 


O-1’  >*2;  (fe-  f  /  S- 

n-o  i 


(IV. 6. 34) 


Equations  (IV. 6. 31)  and  (IV. 6. 32)  are  used  to  obtain  the  gasdynamic  parameters 

at  the  matching  point  X*,  and  thus  the  coefficients  3.  .  f  ,  £  ,  L  and 

0^3  ” 
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Having  computed  all  variables  at  y  =  •  ,  one  can  proceed  to  obtain  the  solu¬ 

tion  at  y  -  t/.  ^  by  applying  the  same  procedure.  The  procedure  is  then  re¬ 
peated  to  cover  the  whole  y  domain. 

t 

iv)  Results  and  Conclusions 

In  order  to  obtain  a  numerical  solution,  it  has  been  considered  that  the  deto¬ 
nating  medium  is  a  stoichiometric  methane-air  mixture,  for  which  ^  =  20.0231, 

"6  3  1.3  and  =  0.033786  (see  section  (IV. 4)). 

The  values  of  the  parameters  of  transport  properties  which  have  been  used  to 
obtain  the  inner  region  solution  are  Re  =  10^,  Pr  =1  and  No  =  1. 

Introducing  the  effect  of  transport  properties  inside  the  inner  region,  the 
density  at  the  center  of  symmetry,  as  expected,  took  the  non-zero  values  shown 
in  Figs.  IV. 6. 3  and  the  temperature  took  finite  values  shown  in  Figs.  IV. 6. 4 
for  different  geometries  and  the  indicated  specific  values  of  R#  ,  Pr-  and 
No  .  The  dotted  range  in  these  curves  expresses  the  continuation  of  the  in- 
viscid  solution  of  the  outer  region  if  transport  effects  were  to  be  neglected. 
Since  the  energy  inside  the  flow  field  is  conserved,  the  solution  of  the  inner 
region  rearranges  the  temperature  distribution  inside  the  inner  region  only  as 
shown  in  Figs. IV. 6. 4.  One  may  expect  that  the  rearrangement  of  such  tempera¬ 
ture  inside  the  whole  of  the  flow  field  can  be  achieved  Dy  considering  the 
transport  effects  in  the  outer  region,  which  would  have  a  rather  small  contri¬ 
bution. 

The  velocity  profiles  for  0.005  £  y  <  0.033,  as  shown  in  Figs.  IV. 6. 5,  extend 
into  negative  values  of  f  indicating  that  the  particles,  after  their  initial 
outward  shift  due  to  the  passage  of  the  front,  returns  towards  the  center  of 
the  flow  field.  This  occurs  due  to  the  relaxation  following  the  strong  explo¬ 
sion  which  causes  nearly  all  the  mass  of  the  blast  wave,  immediately  after  ex¬ 
plosion,  to  be  concentrated  close  to  the  wave  front. 

It  should  be  noted  that  the  pressure  profiles  a're  not  different  from  those  of 
the  inviscid  case.  Figs.  IV. 4.1.  This  result  is  based  on  the  assumption  that 
the  pressure  is  constant  throughout  the  inner  region.  This  assumption  was  up¬ 
held  by  Kamel  et  al. (1977a)  when  he  found  that  the  pressure  distribution  is 


303 


insensitive  to  the  effects  of  both  thermal  conduction  and  radiation. 

The  viscosity  and  conductivity  exponent  <*  is  taken  to  equal  its  self-similar 
values  of  -J,  0  and  1/6  for  planar,  cylindrical  and  spherical  geometries, 
respectively,  whilB  the  radiation  exponent  (3  is  considered  to  equal  7/2,  J 
and  17/6  for  the  same  geometries,  as  given  by  Eq.  (IV. 6. 12).  Therefore,  at  a 
certain  temperature,  the  thermal  conductivity  is  greatest  for  the  spherical 
wave  and  least  for  the  planar  one,  while  the  absorption  coefficient  aS  must 
be  lowest  for  the  spherical  wave  and  highest  for  the  planar  one.  This  means 
that  the  dissipation  effects  are  greatest  for  spherical  waves  and  least  for 
planar  waves. 

Figure  IV. 6. 6  shows  the  thickness  of  the  central  inner  region,  ,  as  a 
function  of  the  parameter  *4  for  plrnar,  cylindrical  and  spherical  geometries. 
As  expected,  the  thickness  of  the  inner  region  has  a  maximum  value  at  =0, 
immediately  after  explosion,  with  highest  heat  dissipation  and  then  decreases 
as  the  wave  front  decays.  At  any  specified  value  of  ,  the  thickness  of 
the  ineer  region  is  greatest  for  the  spherical  case,  with  its  greatest  dissi¬ 
pative  effects,  and  least  for  the  planar  wave. 

Figures  IV. 6. 7  show  the  integral  curves  in  the  f  -  Z  phase-plane,  For  the 
inner  region,  for  planar,  cylindrical  and  spherical  geometries.  As  shown  in 
these  figures,  the  integral  curves,  within  the  range  0.005  <ij^;  0.033,  extend 
into  negative  values  of  f  ,  the  same  tendency  in  velocity  profiles. 

Figures  from  IV. 6. 8  to  IV. 6. 12  give  the  relations  between  the  coefficients 
30  ,  fve  ,  h.*,  »  Pj  and  4  antl  the  Parameter  3  >  respectively,  for  dif¬ 

ferent  blast  wave  geometries. 

Another  important  result  of  this  analysis  is  that,  when  the  mass  integtal,  0^  , 
and  the  energy  integral,  <J3  ,  were  evaluated  with  transport  phenomena  taken 

into  account,  they  differ  very  slightly  from  those  of  the  inviscid  medium. 

This  directly  indicates  from  Eq.  (IV. 5. 14)  that  the  shock  wave  similarity 
front  crajectory  remains  nearly  unaffected  by  transport  phenomena  and  dues 
not  differ  from  that  of  the  inviscid  medium. 

The  sharp  transition  between  the  inner'  legion,  where  tran-opert  pnencmena  are 
most  important,  and  the  outer  inviscid  region  is  as-..cc,iat jo ,  uf  couise,  with  a 


physically  inadmissible  idealization.  However,  smooth  transitions  may  be 
achieved  if  series  expansions  of  the  gasdynamic  parameters  of  higher  powers 
for  the  inner  region  solution  are  used  and  if  the  outer  region  contains  al¬ 
so  the  transport  effects. 


305 


Figure  Captions 

Fig.  IV. 6.1.  A  schematic  diagram  shows  how  the  correct  value  of  the  matching 
point  is  obtained. 

Fig.  IV. 6. 2  A  schematic  diagram  shows  how  the  particle  velocity  at  the  mat¬ 
ching  point  f  is  obtained. 

Fig.  IV. 6. 3  (a),  (b)  &  (c) 

Non-dimensional  density  profiles  of  blast  waves  in  a  detonating 
methane-air  mixture  with  transport  phenomena  taken  into  account 
for  R«  =  104,  Pr  -  1,  No  =  1  and  *  1.3  while  j  “0, 

1  and  2,  respectively. 

Fig.  IV. 6. 4  (a),  (b)  &  (c) 

Non-dimensional  temperature  profiles  of  blast  waves  in  a  detona¬ 
ting  methane-air  mixture  with  transport  phenomena  taken  into  ac¬ 
count  for  R*  =  IQ4,  ?y  =  1 ,  No  =  1  and  =  1.3  while 
$  -  0.  1  and  2,  respectively. 

i 

Fig.  IV. 6. 5  (a),  (b)  S  (c) 

Non-dimensional  velocity  profiles  of  blast  waves  in  a  detonating 
methane-air  mixture  with  transport  phenomena  taken  into  account 
for  Re  =  10^,  Pr  =1,  No  =  1  and  £  =  1.3  while  i  =  0,  1  and 
2,  respectively. 

Fig.  IV. 6. 6  Nar-dimensional  thickness  of  the  inner  region  ,  as  a  function 
of  y  .of  blast  waves  in  a  detonating  methane-air  mixture  with 
transport  phenomena  taken  into  account  for  Re  =  IR^,  Sr  =  1. 

No  -  1  and  If  =  1.3  while  J  =  0,  1  and  7. 

Fig.  IV. 6. 7  i.a),  (b)  &  (c) 

Integral  curves  in  the  F-Z  phase  plane  for  the  inner  region,  for 
different  values  of  ^  .of  blast  waves  in  a  detonating  met  hane- 
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air  mixture  with  transport  phenomena  taken  Into  account  for 
R*  =  10^  >  Pr  =  1,  hfo  =  1  and  K  *1.3  while  j  *»  0,  1  and 
2,  respectively. 

Fig.  IV. 6. 8  The  coefficient  90  as  a  function  of  the  parameter  if  for 
J  «  Dj  1  and  2. 

Fig.  IV, 6. S  Th*  coefficient  lx0  as  a  function  of  Ine  parameter  Jj  for 
j  *  0,  1  and  2. 

The  coefficient  V*  as  a  function  of  the  parameter  y  for 
J  a  C,  'i  and  2. 

I 

The  roefficient  fj  as  a  function  of  the  parameter  y  foi' 

4 

J  *  0,  1  and  2. 

Fig.  IV. 6. 12  The  coefficient  ^  as  a  function  of  the  parameter  y  for 
J  «  0,  1  and  2. 
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